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THE BOUNDARY LAYER DEVELOPMENT IN OPEN CHANNELS 


J. W. Delleur,! J. M. ASCE 
(Proc. Paper 1138) 


SYNOPSIS 


The development of a turbulent boundary layer is analyzed theoretically 
for the case of a steady flow in a horizontal open channel. The problem is 
solved for a rectangular cross-section, but the same scheme could be used 
for other sections. A simplified solution is given for a channel of large 
width. The theory is subsequently compared with experimental results. 


NOTATION 


The letter symbols are defined where they first appear, and are listed in 
the Appendix for reference. 


INTRODUCTION 


The theory is based on the principles of continuity, of constancy of energy 
in the potential flow, and of conservation of momentum in the boundary layer. 
Blasius’ formula is used for computing the friction factor, and the 1/7 power 
law to describe the velocity field. 

An interaction of the boundary layer and the main flow outside the layer is 
demonstrated, with the result that the layer develops more slowly than for a 
flat plate in an infinite fluid under similar flow conditions. For flows close 
to critical the boundary layer may reach a maximum and then decrease. For 
subcritical flow the boundary layer may reach the surface, after which the 
flow is fully developed and remains uniform. 

The experiments have shown that the average behavior of the boundary 
layer is in close agreement with the theory. The tests also revealed 


Note: Discussion open until June 1, 1957. Paper 1138 is part of the copyrighted 
Journal of the Engineering Mechanics Division of the American Society of Civil 
Engineers, Vol. 83, No. EM 1, January, 1957. 


1. Asst. Prof. of Hydr. Eng., School of Civ. Eng., Purdue Univ., Lafayette, 
Ind. 
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the existence of secondary motions which affect the distribution of velocities 
as a result of which the displacement thickness is not constant across the 
channel at a given cross section. 

For flows well in the tranquil range, the theory predicts where the bound- 
ary layer reaches a complete development, which is often of importance in 
certain laboratory investigations. For flows close to critical, the theory 
contributes to a better understanding of the mechanism of the establishment 
of flow over broad crested weirs, for example, and yields a simple way of 
computing a discharge coefficient. 


Background of the Problem 


Some free-surface problems are customarily solved by assuming an ideal 
non viscous fluid. The usual equation for the discharge over 2 broad crested 
weir, for example, is derived in this manner. The neglect of the friction 
may, however, be corrected by considering a layer of shear fiow along the 
boundary. The flow is thus separated into two regions: a boundary layer, 
where the viscous effects are confined, and a region of potential flow outside 
the layer, where the flow behaves as if it were non viscous. The two regions 
of flow are treated separately, as originally suggested by Prandt!* in his 
theory of the boundary layer. 

Consider a steady flow in a horizontal open channel. The flow outside the 
boundary layer is essentially frictionless. Berroulli’s tneorem applies, and 
the specific energy, Ez, ci the potential flow remains constant, that is: 


E. = n+ uj /24 = constant (1) 


where h is the depth of flow, and u, is the velocity outside the boundary layer. 
Furthermore, the principle of continuity, or constancy of the discharge, must 
apply along the channel. 

The boundary layer thickness is designated by §. Two parameters are 
used to describe the boundary layer: the displacement thickness and the 
momentum thickness. The displacement thickness Sy is the distarce by 
which the potential flow has been displaced from the boundary due 'to the 
presence of the boundary layer, and is defined by 


h 
(1- #ydy (2) 


where u is the velocity at the distance y from the boundary. The momentum 
hickness @, is a similar length parameter associated with the momentum 
defect suffered by the fluid because of friction; it is defined by 


(1-4) 4y (3) 


2. Prandtl, L. “Uber Flussigkeiten bei sehr kleiner Reibung” Verhandlung, 
Third International Matmematics Congress, Heidelberg, 1904. 
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The determination of a turbulent boundary layer development by semi- 
empirical methods usually involves the use of the so-called boundary layer 
momentum equation. It is a differential equation derived from the equation 
of motion or making use of the momentum principle.” It was originally de- 
rived by von Karman? for two-dimensional laminar boundary layers: 


where T, is the boundary shear, p is the density of the fluid, and c isa 
friction factor. Assuming no transverse flow, this equation holds for a flow 
in a straight open channel due to the symmetry of the flow. Furthermore, it 
also holds for turbulent boundary layers, by taking u, as the temporal mean 
velocity outside the boundary layer. Expressions for the corrective terms 
ani into account the anisotropy of the turbulence have been derived by 
Ross” and Robertson®, These are shown to be negligible for zero pressure 
gradient. The use of Eq. (4) to describe the conservation of momentum in a 
turbulent boundary layer in open channel flow is, therefore, justified. 
Furthermore, in this particular problem the details of the turbulent flow are 
not known and not required; only the development of the boundary layer is 
sought. 

In order to integrate von Karman’s momentum equation, some relations 
must be assumed or derived between dy, 6 and c. It is known that for 


flows with smali pressure gradients, the shape parameter, k, of the boundary 


layer remains approximately constant, that is 


QO = k = constant (5) 


The only remaining difficulty in applying Eq. (4) is the determination of a 
Suitable expression for the wall shear stress coefficient. For the range of 
Reynolds numbers expected in open channel flow, Blasius formula’ is the 
most convenient. Blasius formula was originally derived for turbulent flow 


in smooth pipes for Re = a € 100,000 as 


f = 0.3164 (6) 


3. See, for example, “Modern Development in Fluid Dynamics,” edited by 
S. Goldstein, Oxford, 1938, p. 131. 

4. von Karman, Th., “Uber laminare und turbulente Reibung,” Zeitschrift fur 
angewandte Mathematik und Mechanik, Vol. 1, p. 235-252, 1921. 

5. Ross, D., “Integration of the Reynolds Equation for Incompressible Turbu- 
lent Boundary Layers,” Ordnance Research Laboratory, The Pennsylvania 
State College, Rept. No. 2, Project NR 062-139; 1953. 

6. Ross, D., and Robertson, J. M., “An Empirical Method for Calculation of 
the Growth of a Turbulent Boundary Layer,” Jour. of Aero. Sc., Vol. 21, 
p. 355, 1954. 

7. Blasius, H., Forschungsheft 131 des Vereins deutscher Ingen. 1911. 
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a where d is the pipe diameter and V the kinematic viscosity, and f is the 
ae pipe resistance coefficient in 


Blasius’ formula (6) may be adapted for use in flows along a flat plate, in 
which case it becomes 


= 0.0225 


§ 


Blasius’ formula is often used in conjunction with the 1/7th power law 
which describes the velocity field in a single function 


for which 


ot (9) 


& 


Blasius’ friction coefficient formula and the 1/7th power law are known to 


be reliable for Reynolds numbers, Re, = 


| s The Concept of Virtual Channel 
| 
v Ee The principles of constant discharge may be combined with the boundary 
{ & layer displacement thickness. Calling q the discharge per unit width, 
. h h 
2 ud u -u d +u h 
2 4 or, by virtue of the definition of the displacement thickness 


q =u, (h -8«) (10) 


- . Eq. (10) states that the actual discharge in a two-dimensional flow is 
8. Schlichting, H., “Lecture Series ‘Boundary Layer Theory,’ Part Il, Tur- 
a bulent Flows,” NACA TM No. 1218, 1949, p. 21. 


it 

L d 2 =f 

where m = a = 

13 
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n equivalent to that of an ideal fluid flowing in a virtual channel of depth h a ‘é 
F This conclusion may be generalized to a three-dimensional channel for which 
ft the actual discharge is that of an ideal fluid in a virtual channel, the bound- 
aries of which are moved inward by an amount oy . The expression of the 
discharge is thus 


Q-Au, (11) 


. where A is the cross sectional area of the virtual channel. i 
Close to the corner there is a zone of mutual interference of the layers ; 
developing along the wall and along the corner. The shape of the boundar a 
layer near a corner has been investigated by Loitsianskii and Bolshakov. . 4 
Their results for a 90° corner and a turbulent boundary layer have been used 
to sketch the edge of the displacement thickness in Fig. 1. The cross section 


of the virtual <hannel is 


A ebh-bS, -2hSy (12) 


The term 0.5 5 — is very small compared to the other terms of Eq. (12). 
In particular, for the experiments of series “A”, to be described subsequent- a 
ly, 0.5 8,2 < 0.001A. Eq. (12) is then simplified to a J 


Aabh-bd, -2hé, (13) 


The Interaction of the Boundary Layer and the Main Flow 


Introducing the expression for the area of the virtual channel (13) into the Lis 
equation of continuity (11), the latter becomes Dos 


Q = ui (bh -2hde) (14) 


The constancy of energy of the potential flow, outside the boundary layer re- : aa of 
quires 


2 


where H is the critical depth. | y 
Eleminating h between (14) and (15), solving for «i and letting q = Q/b, ia 
it follows that 


(16) 


or 


9. Loitsianskii, L. G. and Bolshakov, V. P.,“On the Motion of Fluids in Bound- 
ary Layer Near Line of Intersection of Two Planes,” NACA TM 1308, 1951. I 
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q 
U (17) 


p+ 2h 
b 


Eq. (16) yields the average displacement thickness at any section when the 
specific energy, the surface velocity, and discharge are known, or when the 
surface velocity, depth and discharge are known. This shows the interaction 
of the boundary layer and the main flow, an important feature of the internal 
mechanism of the flow in open channels. It is this aspect that makes this 
boundary layer development different from that along flat plates or airfoils 
in an unconfined fluid. In an open channel the displacement thickness has a 
larger effect on the outside flow because it acts on a finite depth of fluid, 
whereas in an infinite flow, the effect of the displacement thickness is 
negligible as far as the change in the potential flow pattern is concerned. 
Regarding Eq. (17) the difference between the terms of the numerator is 
small, and therefore it yields good results only when they are known ac- 
ecurately. 

Usually the precise measurement of h is difficuit. Eq. (16) may be used 
where the specific energy can be measured with greater accuracy than the 
depth at the station, as is often the case in laboratory work. 


Equation of the Boundary Layer Development 


The equation for the development of the houndary layer is obtained by in- 
troducing the expression for the displacement thickness, a , found above, 
(16), and its derivative into the boundary layer momentum equation (4). The 
velocity u, is taken as the independent variable, and _the ratio of the mo- 
mentum thickness @ to the displacement thickness a is taken as a con- 
stant k, as stated in (5). 

~ It is convenient, however, to use dimensionless parameters. The critical 
depth H = 2 E, is used as the reference length. The following dimension- 
less parameters are introduced: 


bHYgH 


n, € and A are dimensionless displacement thickness, length and channel 
width, respectively; r is a Froude number in terms of the velocity of the 
main flow, and the critical depth; it will be the independent variable. Since 
the criticai velocity is u, = )V gH, it follows that r=1 corresponds to cri- 
tical flow, r < 1 to subcritical flow, and r > 1 to supercritical conditions. 
Similarly, o has the form of a Froude number and has the value of unity for 
critical conditions. The parameter r may be interpreted as the Froude num- 
ber of the potential or main flow, and @ as the Froude number of the total 
flow. 

Introducing the dimensionless parameters of (18) in the equation for the 
displacement thickness (16) and taking the derivative, it follows that 
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where the second form in each case is written to simplify the computation of 
the limit as A approaches  , for the case of a two-dimensional flow, which 
will be considered in a subsequent section— 

The Karman momentum equation (4), in terms of the new variables be- 
comes 


The differential equation governing the growth of the boundary layer is now 
obtained by replacing the values of n and d7 in (21) by (19) and (20): 


_ -3r +r?) } de 


(A +3 -r2)2 


If Blasius’ friction factor (7) is used, it may be expressed in terms of the 


dimensionless parameters (18) as 
jes H 1-3 


which may be introduced in (22) to yield 


i (aks) § -( - 


dr 


(24) 
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Equations (19) and (24) completely define the growth of the boundary layer in 
terms of the discharge parameter g, and the aspect ratio A. The rate of 
growth of the boundary layer may be represented graphically by plotting 
=> 


H against 7, where ~ is a dimensionless length along 


H 


the channel, H V 9H /+ is a Reynolds number in terms of the critical 
depth and critical velocity, and 7 is a dimensionless displacement thickness 
(see Fig. 8). Equations (19) and (24) will be referred to as the three-dimen- 
sional equations of the boundary layer growth. 


Boundary Layer Development in a Wide Channel 
The equations of the boundary layer development may be greatly simpli- 


fied in the case of a two-dimensional flow. When the flow is close to critical, 
and when the channel is very wide compared to its depth, the aspect ratio 


» 3 is large. Letting A or b approach infinity, the equations (16), (17), 


(19) reduce to 


(25) 
(26) 


(27) 
Equation (21) of the boundary layer growth becomes, for A > ‘ 
-(k +) dr (28) 


so that the differential equation (24) of the boundary layer development, in 
which Blasius’ friction factor formula has been used to evaluate c, is: 


“o.0225m'4 


Equations (27) and (28) completely define, in a parametric form, the growth 
of the boundary layer, and yield a family of curves of s ( nan 6 vs. 7 


in terms of the discharge parameter oa, (Fig. 2). Comparing equations (27) 
and (29) to (19) and (24), it appears that the bulk of the solution is given by 
the two dimensional equation, and the effect of the side walls and corners 
appears as a correction. The magnitude of this correction permits one to 
evaluate the effect of the sides and corner, (see Fig. 8). 
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General Properties of the Boundary Layer Development 
in an Horizontal Channel 


For simplicity the general properties of the boundary layer development 
will be analyzed in terms of the two-dimensional theory; they may be ex- 
tended to the three-dimensional theory, in which the effect of the sides and 
corners appear as a correction to the two-dimensional theory. 

At the origin of the boundary layer Se = 0, and Eq. (25) yields a relation 


between o andr: 
3-r* 
(30) 


The plot of Eq. (30) is shown as curve I on Fig. 3. Curve I has two branches, 
the left branch for r < 1 and the right branch for r > 1, which correspond 
to subcritical and supercritical conditions respectively. Therefore for each 
value of the discharge parameter go, there are two possible values of the 
velocity parameter, namely, r, and r3, which correspond to subcritical and 
supercritical conditions respectively. 

Under certain conditions, the development of the boundary layer may 
reach a maximum. Assuming that such a maximum exists, at that point 
dn /dr = 0 for which Eq. (27) yields 


(31) 
The plot of Eq. (31) is shown as curve II in Fig. 3, and r, is the value of r 


which satisfies Eq. (31). Using the definitions (18), it follows that the condi- 
tion (31) is equivalent to ui = gq. But since it was established in (10) that 


q =u, (h- 8s) , the condition (31) yields. 


=/9(h-de) (32) 


This indicates that the flow is critical in the virtual channel of depth h - 5. . 
at the point of maximum development of the boundary layer. The boundary 
layer thus grows to a maximum, producing a throat in the virtual channel. 
This throat is the control of the discharge. When the existence of this maxi- 
mum is physically possible, the discharge establishes itself in such a man- 
ner that the flow is critical at the maximum. 

Eq. (27) and (29), which describe the development of the boundary layer 
have a physical significance only when the boundary layer thickness § is 
less than the depth of water h. That is the limit of validity of the theory 


reached when m a os , or by (15), (16) and (18), when 


r(3-r*) 


(33) 


In the limiting case of the boundary layer becoming tangent to the free 
surface, r =1, and the condition for existence of the maximum of the layer 
is 

b-m (34) 

Beyond its maximum, the boundary layer decreased until it reaches a 


point where the tangent is vertical. This occurs when a5. = 0 or by (29) when 
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ali 2k+i +i r*)r (35) 
2 2(k+!) 


The plot of Eq. (35) is shown as curve III in Fig. 3, for k =4 , and where P 


7 r2 is the value of r which satisfies (35). The section of the curve of the 


growth of the boundary layer, between the maximum and the point of vertical 
tangent, that is for r: < r < rz, is not accurate, because of the rapid in- 
crease in velocity between r; and rz which requires a rapid change of sur- 
face slope. This would necessitate taking into account the curvature of the 
flow. This situation may occur when the channel ends with a free overfall. 
Due to the confinement of the flow and the effect of the free surface, the 
rate of growth of the boundary layer in open channels is slower than along a 
flat plate in an infinite fluid. For the purpose of comparison the development 
of the boundary layer along a flat plate is shown in Fig. 2,10 


Laboratory Set-Up 


The general layout of the laboratory equipment is shown in Fig. 4. The 
main feature of the equipment is the test channel which is 15' - 7-1/2" long 
and 13-11/16" wide. It is part of a closed circuit formed by the reservoir, 
the pump, the constant head tank, the entrance section with gravel baffles 

and screens, the test section of the channel, the exit of the channel with a 
measuring weir, and finally, a return channel leading back to the reservoir. 

The measuring equipment consisted of a pitot tube with its traversing 
mechanism, a micromanometer, and point gages. The pitot tube consisted of 
a hypodermic needle inserted in a streamlined lucite stem. Two sizes of 
needles were used: 0,027" O.D. and 0.060" O.D. Almost all data presented 
herein were taken with the 0.060" tube, which had an inside diameter of 
0.041"". The traversing mechanism controlled the vertical position of the 
pitot tube to the nearest 0.001"', and the micromanometer was designed to 
indicate pressure differentials of 0.001" of water. 

The discharge was controlled by maintaining a constant pool elevation at 
the entrance section. The elevation was measured with a point gage which 
permitted the reproduction of the pool elevation for consecutive experiments 
on the same flow. 

The investigation being focussed on the development of a turbulent bound- 
ary layer, it was deemed desirable to ensure that the transition from a 
laminar to a turbulent flow in the boundary layer always occurred at the 
same point. The point selected was the end of the curved inlet, so that the 
boundary layer would be turbulent throughout the length of the channel proper. 
A 1/16" wire was found convenient for the range of experiments. 


Experiments 


Four series of experiments were carried out, with the following charac- 
teristics: 


: s/ ‘i 10. Schlichting, H., “Lecture Series ‘Boundary Layer Theory,’ Part II, Tur- 
i rans bulent Flows,” NACA TM 1218, 1949, Eq. (17.9). | 
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Type of Total Discharge Aspect 
Series Flow Energy Parameter Ratio 
E 


Subcritical .218' .84 7.85 
Subcritical .259" 4.41 
Subcritical .2505' .855 6.83 
Supercritical -6835' .158 2.49 


ASCE 


Series “A” is the most complete of the four series of experiments. It was 
designed to check thoroughly the theory for a particular flow condition. 
Series *B*, “C” and “D® are more qualitative and were intended to explore 
other flow conditions, Velocity profiles were taken at the center of the chan- 
nel (z =0) and at z=+.1'; *+.2"; +3'; +4", where z is measured normal to 
the axis of the channel, at 4 stations: x =52"; 99.75"; 135.75" and 171.5", 
where x is measured from the beginning of the channel. In addition, center- 
line velocity profiles were taken at x = 4"; 25"; 63.875"; 81.75"; 123.5" and 
153.69"". Typical data are presented in Figs. 5 to 7. 

The displacement thicknesses Se were computed by numerical integra- 
tion from the velocity profiles. Where several profiles were taken at the 
same station, the average Sa was computed (a weight of 1 was given to each 
value of be , except those pertaining to outer profiles which were assigned a 
relative weight of 1/2). These average displacement thicknesses were then 
expressed in dimensionless form and compared with the values predicted by 
the theory. The displacement thicknesses were also evaluated from Eq, (16). 
The results are summarized in Table 1, and in Fig. 8. 

Series *B” has a low value of go, indicating a very tranquil flow, the ve- 
locity being far from critical. The boundary iayer in this case would pierce 
the surface before reaching its maximum, The experimental channel was too 
short to reach the point of full development of the flow. (See Fig. 9.) 

Series “C” was intended for the check of formula (16). Thé local depths 
were measured with a point gage indicating the nearest 0.001". The dis- 
placement thicknesses, computed by this formula, are compared with the ex- 
perimental data, and the theoretical curve of the boundary layer growth in 
Fig. 9. 

Series *D” was intended to check the theory in the case of supercritical 
flow. The supercritical condition was obtained by placing a sluice gate at the 
entrance of the channel, The rough surface for this type of flow made the ob- 
servations difficult; the results are shown in Fig. 9. 


Analysis of Data 


Examinations of Fig. 8 reveal the main conclusions to be drawn from the 
experimental investigation. The average displacement thickness confirms 
the three-dimensional theory. The shape of the theoretical curve of the 
boundary layer development, the existence of the maximum point and the sub- 
sequent decrease of the boundary layer thickness are also confirmed, The 
curve through the average experimental displacement thickness is consistent- 
ly slightly below the theoretical curve. The individual experimental data from 
points across the channel at a given station are scattered about their mean, 
This scatter is attributed to the presence of secondary motions, The theory 
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thus predicts the mean behavior of the boundary layer, but does not show the 
detail of the boundary layer distribution as influenced by the secondary flow. 
Motion pictures of tracer particles have shown the existence of these second- 
ary flows. 

The boundary layer experiments also showed that formula (16), which also 
evaluates the mean behavior of the boundary layer, yields results which are 
in close agreement with the theory. 

The experiments show that the 1/ 7th power law describes adequately the 
average velocity field. The 1/7th power law has been generalized to a power 


n 
formula of the type 4 = (2) , for which 


The experimental data indicate that n varies from a maximum of = to 
1 
a minimum of 770° the average n for the experiment was 685° The nu- 


merical computation could have beer carried out with this value of n instead 
of 1/7. However, this slight deviation does not affect the results, because in 
the equation for the growth of the boundary layer n occurs indirectly in m¢+ 
and k, which are not very sensitive to small changes in n. 

It may also be remarked that the equations (19) and (24) governing the de- 
velopment of the boundary layer do not contain explicitly the length of the 
channel. This is due to the fact that the specific energy Eg, the discharge 
parameter g and the length L of the channel are not independent, For ex- 
ample, if the energy Eg and the length L are specified, the discharge para- 
meter go is determined. In other words, o is not a free parameter when 
E, and L are given. This fact may be used to analyze, for example, how 
the discharge establishes itself over a broad crested weir, and yields a 
theoretical way of computing a discharge coefficient, thus relating the fixed 
geometry of this hydraulic device and the discharge when the boundary layer 
is the governing factor.16 The governing influence of the boundary layer on 
the rate of flow over embankment shaped weirs has also been substantiated 
in the comprehensive experiments carried out by Sigurdsson,!7 

He also confirms the use of the 1/7th power law to describe the velocity 
distribution, It must be pointed out that the application of the boundary layer 
theory to broad crested weirs has certain limitations. In fact, other phe- 
nomena may play an important role in the hydraulics of broad crested weirs, 
such as the curvature and possible undulation of the surface. There may 
also be a separation at the entrance if it is not rounded, which would have a 
predominant effect. The theory is applicable only where the flow is essen- 
tially parallel and free from undulation and separation. 

For flows well in the tranquil range, the boundary layer grows until it 
reaches the surface, after which the flow becomes uniform. The determina- 
tion of the distance to the point of uniform flow in intake section of open 


16. Delleur, J. W., “The Boundary Layer Development on Broad Crested 
Weirs,” Fourth Midwestern Conference on Fluid Mechanics, Purdue Uni- 
versity, 1955. 

17. Sigurdsson, G., “Discharge Characteristics of an Embankment-shaped 
Weir,” M.S. Thesis, Georgia Inst. of Tech., June, 1955. 
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channels is important in certain laboratory investigations, 18 Or, inversely, 
the discharge which will yield a uniform flow in an open channel may be re- 
quired. 


CONCLUSIONS 


The development of a turbulent boundary layer is treated for steady flow 
in a rectangular open channel. A simplified solution is given for a channel of 
large width. 

The theory shows an interaction of the boundary layer and the main flow 
outside the layer, with the result that the layer develops more slowly than for 
a flat plate in an infinite fluid under similar flow conditions. For flows close 
to critical the boundary layer may reach a maximum and then decrease. For 
flows less than critical the theory predicts where the boundary layer reaches 
the surface. 

. The experiments have shown that the average behavior of the boundary 
layer is in close agreement with the theory. The tests revealed that there 
are secondary motions which affect the distribution of velocities. Asa re- 
sult the displacement thickness is not constant at any one cross section, but 
the average displacement thickness is in close agreement with the theory. 

Two applications of the theory are mentioned: the distance to uniform 
flow in the intake section of open channel, which is sometimes an important 
factor in laboratory investigations; and the computation of a theoretical dis- 
charge coefficient for flows over horizontal broad crested weirs. 

From a broader point of view, the consideration of the boundary layer in 


open channel opens new horizonts for a better understanding of the internal 
mechanism of this type of flow. 
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o7 97 


n 
v 
§ 

p 


friction coefficient for pipe flow 
depth of flow 
critical depth 


: shape parameter of velocity profile 


: length of zone of interference of bottom and wall boundary layers 


ratio 
discharge per unit width 
total discharge 


: velocity parameter, dimensionless 


Reynolds number 


: longitudinal velocity 


velocity at the edge of the boundary layer 


: distance along the channel 


boundary layer thickness 
displacement thickness 
ratio da/t 

ratio b/H 

kinematic viscosity 
ratio x/H 


: density 


discharge parameter, dimensionless 


boundary shear stress 
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FIG. T. VELOCITY MEASUREMENTS AT X =99.75" 
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APPROXIMATION TO PLASTIC BEHAVIOR OF CIRCULAR MEMBRANES 


N. A. Weil,! J. M. ASCE 
(Proc. Paper 1139) 


SYNOPSIS 


An approximate theory is presented for the plastic behavior of clamped 
circular membranes subjected to lateral hydrostatic pressures. It is as- 
sumed that the actual stress-strain curve of the material can be well repre- 
sented by a parabolic expression. The von Mises-Hencky relations are taken 
to govern the plastic behavior of the material; suitable approximations are 
then adopted for the biaxial stress-strain and strain-displacement laws, The 
derivation is based on infinite power series expansions, wherein terms con- 
taining the higher powers of small quantities are neglected. 

The predictions of this approximate theory are compared with experiment- 
al results obtained from bulge tests of annealed copper diaphragms. Theory 
and experiments are in satisfactory agreement in the low and medium pres- 
sure ranges. In the vicinity of instability the theory gradually deviates from 
experimental observations; this is accounted for by the increased contribution 
of the omitted terms near the bursting pressure of the diaphragm. 


NOMENCLATURE 


Radius of diaphragm, in = h/h = Dimensionless thinning 
Thickness of diaphragm, in r = Radial distance from pole, in 
True stress, psi p =r/a = Dimensionless radial 
Natural strain distance 
Octahedral shear stress, psi w = Vertical deflection, in 
= Octahedral shear strain w = w/a = Dimensionless deflection 
Hydrostatic pressure, psi u = Horizontal displacement, in 
= Pa/h = Relative pressure, psi v =u/a = Dimensionless horiz. 
displacement 
m = Strain hardening exponent 
A = Hardness factor, psi 


Note: Discussion open until June 1, 1957. Paper 1139 is part of the copyrighted Jour- 
nal of the Engineering Mechanics Division of the American Society of Civil Engineers, 
Vol. 83, No. EM 1, January, 1957. 

1. Design Engr., M. W. Kellogg Co., New York, N. Y. 


1139-1 


[2 
Paper 1139 EM 1 
° 
; 
4 
«4 
4 
— 
“| 


1139-2 EM 1 January, 1957 


Subscripts r, 9, z denote entities acting in meridional, circumferential 
and normal (i.e. principal) directions. Bars above symbols represent enti- 
ties in the original (undeformed) plate; the subscript o refers to variables at 
the center (pole) of the bulge. 


INTRODUCTION 


The problem of the plastic behavior of clamped membranes is important 
from two structural viewpoints: to determine the load-deflection relations in 
the plastic domain and to obtain an estimate on the ultimate lo-d carrying 
ability (or the bursting pressure in the case of fangible disks) of circular 
plates. Rheologically this problem is used to examine the validity of various 
plasticity laws over a wide range of strains. 

All of the aforementioned aspects of this problem have received the atten- 
tion of investigators. An analytical treatment, based on the Hencky- Mises 
plasticity laws and strain-displacement relations valid for large deformation 
in the plastic domain, was recently presented in the literature.(1) This solu- 
tion, hereinafter referred to as the “exact theory,” is in good agreement with 
experimental results, and predicts the significant values at instability with a 
fair degree of accuracy. The solution itself, however, is cumbersome, re- 
quiring a numerical iteration scheme between two simultaneous integral equa- 
tions. 

To avoid the difficulties inherent in an exact solution, Hill (2) proposed an 
approximate solution, whose results are given in closed form. This solution 
is based on the Levy- Mises plasticity laws (i.e. strain-rates rather than 
strains are considered), and is valid only for materials whose stress-strain 
curve can be well approximated by a linear relationship. 

Most materials, especially annealed metals, do not possess a constant rate 
of strain hardening. For materials with a decreasing rate of strain harden- 
ing the Hill-theory, as will be shown later, may be insufficiently accurate to 
be used even as a crude approximation. In this case the assumption of a 
power-function (parabola) for the stress-strain relation would lead toa 
closer correlation with the physical behavior of the structure investigated. 

This assumption forms the core of the current work. To comply with the 
nature of the problem, strain-displacement relations valid for large deforma- 
tions were adopted, although their form was somewhat simplified from that 
used in the exact solution (1) to lessen the complexity of the derivation. The 
material is also assumed to be insensitive to rate of straining. The laws of 
static equilibrium are observed and the plasticity laws are approximately 
satisfied. The solution is obtained in the form of successive approximations, 
in which higher powers of small quantities are ignored. The current work is 
restricted to the full development of the second approximation. The dimen- 
sionless notation used in the derivation allows for an immediate application 
of the results to diaphragms of various radii and thicknesses. 


Solution of the Approximate Theory 


The von Mises flow condition postulates that plastic deformation will occur 
when the octahedral shear stress 


- 0, 


| 
i. 
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exceeds a given limiting value. The magnitude of flow is measured by the 
octahedral shear strain 


Y2 (2) 


2 2 2 
2 - + (e,- + - 


where the yield curve 7 = 7 (Y) is assumed to be known. Assuming volume 
constancy 


(3) 


the flow condition for a uniaxial stress 01 and a longitudinal strain € 4 be- 
comes 


In the case of biaxial stress (0; = 9, €4 = €9), the flow condition modi- 
fies to 


(5) 


Eqs. 4 and 5 show that plastic deformation is governed by identical rela- 
tions under uniaxial and biaxial stresses, if in the latter case the negative 
through-thickness strain is taken as a measure of the amount of straining. 

Assume now that uniaxial stress-strain relation for a given material fol- 
lows the power-function law 


2(i-m)/2 a 73 


where m is the strain hardening exponent and A is the stress value (hardners 
factor) at = 1. 

In a clamped diaphragm subjected to lateral pressure, the stress distribu- 
tion varies from a balanced biaxiality at the crown to a biaxiality ratio of 1/2 
at the edge. It is, therefore, a reasonably good simplification to assume that 
the yield criterion will be given by 


= A(-e,)” (7) 


4 

+ = O 

(4) 
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Let the finite strain-displacement relations be given by the approximate 
dimensionless form: 


=%/p (8) 


where primes denote differentiation with respect to Pp. 

Hill’s first approximation, (2) postulating a constant stress throughout the 
diaphragm, shows that the through-thickness strain will be given by (to this 
approximation) 


(1 - p~) (9) 


Adopting this criterion as the starting point for the second approximation, 
Eq. 7 assumes the form 


The equilibrium of normal forces yields [see Eq. 13 of Ref. (1)] 


ul 11 
% s (2-3 )>, * por ( ) 


One can eliminate 0g from Eqs. 10 and 11 to obtain 


2 1+ (1 - 12 
Furthermore, Eq. 9 and the last of Eqs. 8 can be combined to give 
(13) 
1/1 
which substituted into Eq. 12 leads to 
2\m 
1 2 lim, 2m (1 - 1 


Eq. 14 is a first order differential equation with variable coefficients, p 


a 
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being the independent and ¢, the dependent variable. The explicit solution of 
this type of differential equation is well known; however, the integrals appear- 
ing in this solution cannot be evaluated in closed form. 

Therefore, in what follows, infinite series expansions will be used, unless 
the problem lends itself to a solution in closed form. In doing so, only those 
terms of the infinite series are retained which involve powers of the variables 
no higher than p°, wo° or m2. 

These assumptions merit more detailed consideration. The range of the 
dimensionless radial distance is 0< Pp < 1. In the vicinity of the clamping 
edge, p < 1, stresses and strains are relatively small, so that omission of 
terms containing the sixth or higher powers of P should not introduce serious 
errors. 

For annealed materials one commonly finds that m < 0.4 and for severely 
cold-rolled materials m < 0.1. Hence the error due to the omission of terms 
in m3 or higher will generally stay below 7 per cent; in most instances the 
error will be much less, remembering that most series expansions have in- 
verse factorial coefficients. 

Wo varies from 0 at the beginning of the experiment to a maximum possible 
value of 1; most engineering materials fail at about Wo < 0.6. At low pres- 
sures the contribution of terms of order wo or higher is entirely negligible. 
At increasing pressures, however, the error due to the omission of these 
terms becomes more significant. This will lead to a gradual deviation of the 
actual values from the results predicted by this approximate theory, the er- 
ror being largest at instability. 

Proceeding on the principles set down above, and disallowing infinite 
stresses at the pole, the solution of Eq. 14 yields o,;; 9g can then be obtained 
from Eq. 10. The stresses will therefore be given by the expressions 


2, 2 2. 4 
° - (0. 5m+w,, - [ 0. 5m(1-m) mw, 


% 


2 4 
= ™ {2 - - [ 2. 5m(1-m) +m, + ] 


The relationship between pressure and vertical deflection can be obtained 
in the following manner. Eq. 15 of Ref. 1 can be modified to 


and expanding 


(16) 


The term containing p 5 has been omitted in Eq. 16 because otherwise the 
solution would become extremely complex. This step, however, is justified, 
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since the numerical coefficient of this term is exceedingly small (1/128). 
Eq. 13 can be readily solved for n subject to the boundary condition that 
n = 1at p=1. Expanding this solution in series one obtains 


2,7 + au (17) 


Inserting now co , from the first of Eqs. 15 and 7 from Eq. 17 into Eq. 16, 
results, to this approximation 


+ [ + mwo* + (2-200? bag | 
6 


+ + + [ - (1-3m)w,* - (9m | 
Pp 3 
[ (0.5+m) = 2 + 4] 
° 2m 
° 


For a constant pressure the only independent variable in this equation is 
P. Hence it can be integrated to yield 


4 


6 
+ [ m(2+0. 5m) + + | 
Aw, 


- + (-3m)u,” - (6-9m)u,"|E + 


8 3 
2 Aw, 


In Eq. 19 wo appears as a parameter. Its value can be obtained from the 
boundary conditions: 
W =Wo at 
w =0 at p 
Substitution of the boundary values and elimination of C from the resulting 
two expressions of Eq. 19 leads to 
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[2 + 1, 25m*+ O.5m + 5m+2m” * + (8+10, 5m | 
2 


+ [ 1+ 0.3056m + 0, 0278m* + (1. 5+0, 5m+0, 0556m")w,* 


+ (1, 1111+0, 3889m+0, 5m dw P -w, = 
° 2m, 2m ° 
2 dw, 


The above expression is cubic in the relative pressure, p. Treating this 
quantity as the independent variable (w, being fixed for a given pressure) and 
using the Newton-Raphson method,(3) Eq. 20 can be solved to yield 


p= [2 - 0.3056m + 0, 0656m* - (2, 5-0,3889m-O, 118m” + (21) 


+ (6,139+0, 5332m2. | ra 


Having obtained the value of C from the boundary conditions and p from 
Eq. 21, one can return to Eq. 19 to obtain an explicit expression for the 
dimensionless deflection: 


w [ 0,3056m 0,0656m~ + (0,5+0,222m-0, + 


(1, 8611-0, 69982, 5180 | p* + 


+ [ 0.05560 0,0108m* + (1. 1111-0, 63731" + 


Furthermore, Eq. 21 can be inverted to solve for w, in terms of p. There 
results 


1 + 0,3056m = 0,5833m~ + (2,52,333 ) 
(23) 


2+m 


+ (12, 611-33, 99m+79. 56n°) ( 
2 


Adopting now the Hencky deformation theory for the yield criterion 


3 at 
+ 
4 
a) 
(22) 
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(24) 
i,j,k=1,2,3 iA 


and using the law of volume constancy, Eq. 3, the following identity can be de- 
rived to relate principal strains and stresses 


25. - 9 
(25) 
~ % 
and by substitution of Eqs. 15 and the second of Eqs. 8 
[2 + Gm-bu,7)p* + | 


The first of Eqs. 8, in turn, can be rewritten as 


which, upon substitution of Eqs. 22 and 26, leads to the first order differential 
equation 


2 
- E + Gm-bw,”)o + ] 2. 


1 - 0,6111m + 0,22,,5m* - (1+0, 1389m-0, 569m" - 


(3.472-1, 621 9m+5. 235m” ‘ 


[= ~ « (2-4,.667m+2. 705m") (27) 


| 


[0.3333 + 0,2130m* + (7m9, 222m" )w * 


+ 


- (10, | “| 


Solving Eq. 27 subject to the boundary condition that v = 0 at P = 1 pro- 
vides the expression for v. Then, by the second of Eqs. 8, 
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(9.75-18.079m-113. 65m" ‘ 


(28) 
- (10.778-0. 9992m-130. 4.2m“) | + 
2 2 2 
+ 10,4444m+0,6590m* - (1.6667m-11.134m ) 
+ (0. 7222-28, 44m+103, 12m2) | 
and, from Eq. 26 
{2 - 0,9167m 0,5471m* + (1, * - 
~ 65m" + 
+ [4-290 = (6,5-15.417m-3. 106m~)w * - — 
(29) 


(21+29, 92m-199. osm | 


+ %, 909m" - (29. 33m-69, 2m” + 


The through-thickness strain, € 7, can be derived by simply adding Eqs. 28 
and 29, and providing this sum with a minus sign, according to Eq. 3. The 
polar stress and the in-plane polar strain can be derived from Eqs. 15 and 
28, respectively, by setting Pp = 0: 


o (30) 


° 


6 (31) 
- (9.75-18. 07 9m-113. 
This completes the solution to the second approximation. The next 
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approximation would demand that Eqs. 15 be reintroduced into Eq. 10 in order 
to obtain a more defined relationship for the stresses. 


Comparison of Approximate Theory with Experiments 


The experiments conducted on annealed copper membranes of two thick- 
nesses were described in a previous paper.(1) 

In brief, annealed electrolytic copper plates of two thicknesses were 
tested: two “thin” plates with h = 0.02087 and 0.02020 in., respectively, and 
a “thick” plate with h = 0.1303 in. These plates were clamped into a device 
with a 11.5 in. aperture, rendering the radius a = 5.75 in. 

Pressure was applied in predetermined steps to one side of the plates 
during the tests. After reaching a desired level, the pressure was released 
prior to any measurements, in agreement with theoretical considerations 
which ignored the contribution of elastic deformations. The initially flat plate 
deformed into a bulge as the test progressed; physical measurements (apart 
from pressure values) consisted of strain readings with the aid of a photo- 
grid (20 lines per inch) using a microscope of about X9 magnification as well 
as clip strain-gages with a range of 0.20 in/in total strain capacity. Vertical 
and horizontal displacement were measured by means of a 0.0001" Ames dial 
mounted on a specially constructed bridge. These measurements were car- 
ried out over two mutually perpendicular radii of the plate, with clip gages 
applied to both the top and bottom surfaces. 

In order to bring all tests to a common denominator, test results were 
plotted in terms of the “relative pressure,” p. Conversion from actual pres- 
sure, P, to relative pressure was made through the relationship p = Pa/h. 
Calculated numerical results are also expressed in terms of the relative 
pressure, p. 

The stress-strain curve for this material is shown in Fig. 1. As can be 
seen in this figure, the actual stress-strain curve can be approximated by an 
expression of the type of Eq. 6, if the material constants are taken as 


A = 65,500 psi ; m = 0.375 


Upon substitution of these numerical values into Eqs. 23, 31 and 30 in suc- 
cession, the approximate theory yields the following expressions for the di- 
mensionless crown deflection, polar strain and polar stress for this material: 


° 


2132 17.38 
L282 0.5714 , 52808 51,7143 , 222385 (32) 
10" 10 


+ 1.252%," + (33) 


0.75 (34) 
* BL 


For given magnitudes of the relative pressure p, one can now calculate the 
dimensionless polar deflection from Eq. 32; in possession of this quantity the 
strain and stress at the crown are obtained from Eqs. 33 and 34 respectively. 
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The results of these calculations are shown in Figs. 2, 3 and 4,* along with 
the experimental values ovtained from the tests. The agreement between the 
approximate theory and experiments is quite good at low and medium pres- 
sures. Near instability the agreement becomes less satisfactory; vertical 
deflection and stresses become too low, whereas predicted strains are ex- 
cessive. Furthermore none of the curves exhibit a vertical rise, which 
characterizes the physical behavior of the membrane at instability. 

Some of the reasons for the deviation of this theory from experiments in 
the vicinity of bursting have already been mentioned in the text. Among these 
are the approximate nature of the assumptions incorporated in Eqs. 6, 7, 8 
and 9, the omission of higher order terms from the infinite series expansion 
and the fact that this solution was carried only to a second iterative step. 

The approximation involved in Eq. 6 could be improved by adopting a stress- 
strain expression of the type 


Ae™ + Be™ 


In view of the fact that even the simpler assumption for the stress-strain 
curve is quite close to the actual physical properties (see Fig. 1), this im- 
provement is not felt to be significant. Rather, it seems that the approxima- 
tions involved in Eqs. 7 and 8, and the sizeable contribution of the terms of 
order w & are the chief cause of the discrepancy near the limiting p1essure. 

As mentioned before, a different approximate theory was previousiv pro- 
posed by Hill.(2) It is based on the tenet that the actual stress-strain relation 
can be well approached by a linear expression of the type 


Y(1+He ) (35) 


To the second approximation this assumption yields the following expres- 
sions for the significant variables: 


(6H-10)w,? + lw, = p/t (36) 


(37) 


(38) 


o = *)Y 
° ° 


For a given relative pressure level p, w,_ is obtained from the cubic ex- 
pression given in Eq. 36 by some suitable method, say the Newton-Raphson 
method. € 9 and ¢o can then be readily determined explicitly from Eqs. 37 


and 38. 
For the annealed copper used in this comparison, the material constants 


appearing in Hill’s theory were assumed to be 
Y = 15,000 psi ; H = 6.1 


* In Fig. 2 the dimensionless polar deflection, w 9, has been reconverted into 
actual center deflection, w. 
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The linear stress-strain curve based on these assumptions is also pic- 
tured in Fig. 1. It may be added that this line was obtained by choosing the 
best linear fit to the real o - € curve over the strain range 0< € < 0.32 en- 
countered in the bulge tests. Finally the vertica) center deflection, strain 
and stress predicted by Hill are shown in Figs. 2, 3 and 4 in dashed lines. 
The vertical center deflection and stress, Figs. 2 and 4, in the low and medi- 
um pressure range are predicted as well by Hill as by the theory described 
in this paper. At high pressures the Hill theory is less apt to follow actual 
physical behavior, insofar as these variables are concerned. For the polar 
strain, Fig. 3, Hill’s approximation predicts a maximum, whereafter the 
polar strain is seen to diminish with increasing pressures. This is contrary 
to the physical behavior of disks subjected to lateral hydrostatic pressure. 

The investigation of the instability condition and limiting pressures will be 
dealt with in a separate paper. 
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SUMMARY 


In current design practice it is generally assumed that a truss is equally 
safe against buckling in and out of its plane if its members have the same 
slenderness ratio in both directions. In this paper it is shown that the argu- 
ment may be inapplicable to the customary rigid jointed type of structure, In 
fact it is demonstrated that this kind of truss, designed under present code 
requirements, may exhibit a much smaller factor of safety against buckling 
normal to its plane than in its plane. 

The analysis of the lateral stability proceeds from the establishment of 
the basic equations to the derivation of a determinantal form of stability 
criterion. An alternate approach, pursued further in this paper because of 
its reduced numerical difficulty, utilizes a generalized moment distribution 
technique in the analysis of the response of the truss to applied unit couples 
whose planes are perpendicular to that of the structure, This results in the 
development of a *series criterion” which represents an extension and gen- 
eralization of a similar criterion first proposed by Lundquist. 


INTRODUCTION 


In the design of a compression member of a steel truss, it is customary, 
under present code requirements [1],4 toconsider the allowable axial stress 
a function of the largest slenderness ratio alone, without regard to the condi- 
tion of elastic restraint at the joints, It follows therefore that an “efficient” 
member should display approximately the same slenderness ratio with re- 
spect to both principal directions; this, it is felt, should insure an equal 
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factor of safety against buckling both in and out of the plane of the structure. 

It is fairly obvious that this kind of reasoning cannot be applied to a large 
and common class of trusses without major modification. Even in the event 
of perfect pin connections at the joints, the argument is valid only if the cen- 
troid and shear center coincide. This is often not the case. The customary 
double-angle type of section, for example, buckles laterally in a mode which 
couples bending and twisting, and under an axial load which may be appreci- 
ably smaller than the Euler load, which is based on a bending mode alone, 

If the bars are held rigidly at the joints, this discrepancy is likely to in- 
crease further. When buckling occurs in the elastic region (that is, for suf- 
ficiently slender elements) the critical compressive stress depends to a large 
extent on the degree of end restraint. This restraint, which is imposed on 
one member by the adjoining members, is generally less effective relative to 
lateral buckling that in relation to buckling in the plane of the structure. This 
is especially true of double-angle sections, whose torsional stiffness is negli- 
gible compared with their bending stiffness. It is therefore not implausible 
to expect trusses composed of such members to exhibit a greater tendency to 
buckle normal to their plane than in their own plane. 

It is the purpose of this paper to investigate the conditions under which 
this type of instability will occur. In order to facilitate the analysis, a num- 
ber of assumptions are made which, it is felt, are reasonable without reduc- 
ing the generality of the approach to a significant extent. It is assumed, for 
€aaw.ple, that the members of the truss are joined rigidly at the panel points; 
the effect of the gusset plates is therefore disregarded. It is postulated fur- 
ther that the panel points are held rigidly against lateral translation by means 
of lateral cross bracing, which, however, does not contribute to the rotational 
joint stiffness. Since all members are composed of two angles, the further 
assumption that the shear center itself is held against lateral movement at 
the joints appears not unwarranted in view of the fact that the cross bracing 
is ordinarily connected to the outstanding legs. 

The analysis in the form presented here is not confined to double-angle 
sections; it is restricted, however, to sections whose warping stiffness may 
be considered negligible in comparison with their torsional and bending stiff- 
nesses, An extension to sections exhibiting appreciable warping resistance, 
such as, e.g., I-beams, may be effected at the expense of some of the nu- 
merical simplicity. On the other hand, this does not appear worthwhile since 
members of this type are rarely slender enough to be in danger of elastic 
buckling, at least so far as common structural practice is concerned. 

It is finally assumed that, prior to lateral buckling, all members are sub- 
jected only to axial forces through their centroids. This implies, of course, 
that all external loads are applied at the joints, that all bars are initially 
straight, and that the effect of the joint translations within the plane of the 
framework may be neglected. Actually, “secondary moments” are inevitable 
and are likely to affect the condition of lateral instability to a significant ex- 
tent, Nevertheless, this omission is not considered serious, as will be dis- 
cussed in the Conclusion, 


Establishment of Basic Equations 


In what follows, let a representative bar between joints i and j buckle 
laterally as shown in Fig. 1. The differential equations of equilibrium 
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governing the displacements (u,v) of the shear center S and the rotation B 
are well-known and are given, for example, in [2], p. 132. For the case under 
consideration, in which the warping constant is assumed to vanish and in 
which, by virtue of the symmetry of all sections, buckling in the plane of the 

j framework is uncoupled and may therefore be ignored, they read: 


Ely ulY + Pu" + P Yo B" = O (1) 


= 0 (2) 


+ 


In these equations, the axial force P is considered positive in compres- 
sion, E and G are the usual material constants, while Iy and K are, re- 
spectively, the moment of inertia about the y-axis and the torsional stiffness 
constants,» The radius of gyration about the shear center is designated by 
Pp, and y, is the coordinate of the shear center relative to the centroid C. 
Primes, as usua!, indicate differentiation with respect to the independent 
variable (z) along the axis of the bar. 

Eqs. (1) and (2) can be simplified somewhat. If the dimensionless quantity 


uw is defined by 


(3) 


it follows, from Eq. (2), that 


and, upon substitution in Eq. (1), that 


Ely ul¥+Nu" = 0 (5) 


in which N is an *equivalent” axial force given by 


N = P (1 (6) 


The general solution of these linear homogeneous equations is readily es- 
tablished and involves six constants of integration. Basically, this process 
should now be repeated with each member. There results a large number of 
constants of integration, for which an equal number of linear homogeneous 
boundary conditions is available. The problem is to determine the smallest 
load for which the ensuing system of linear homogeneous equations admits a 
non-trivial solution; in this fashion, a typical determinantal stability criterion 
is established, 

This represents a task of formidable numerical proportions, It can be 
simplified greatly, however, through a systematizing procedure which is 
analogous to the establishment of modified “slope-deflection® equations. To 
this end, the bar rotations 8 and @ at the joints (representing twisting and 


5. To be on the safe side, K for the section is assumed to be twice the value 
of the individual angle. 
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a rotations, respectively) are introduced as independent variables, and 
qs. (4) and (5) are solved subject to the boundary conditions. 
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u(0) u(L) 
u'(0) u'(L) = -53 (7) 
B(O) = B(L) = 


in which L represents the length of the member. The choice of the sign con- 
vention is governed by considerations of consistency for the developments 
below. 

The functions u(z) and 8 (z) are thus established in terms of the four end 
rotation 2; their form is of no direct interest here and has therefore been 
omitted. However, it is now possible to express the twisting moment T and 
bending moment M at the joints in terms of these rotations. This is ac- 
complished by letting 


(GK - Pp®) B'(0) - Py, u'(0) 
(GK - Pp®) B'(L) - Py, u'(L) 


Ely u"(0) 


Ely u"(L) 
In view of these equations, the moment-rotation relationships take the fol- 
lowing simple form in the case of prismatic bars. 


= (Byy - Ce B53) (9) 


with Ty and Mii obtained by inverting the subscripts (i,j) on the right side 
of Eqs. (9). It is of interest to note, and is physically plausible, that the 
choice of the shear center as reference point has led to an uncoupling of 
twisting and bending, so far as the joints are concerned. From this it is not 
to be inferred, however, that there is no coupling all along the bar; on the 
contrary, pure bending at the joints, for example, produces bending and 
twisting, unless the centroid and shear certer coincide. 

In Eqs. (9), the “stiffness factors” S; and S,, and the “carry-over fac- 
tors” C; and C,, are given by the following equations: 


Cy - 1 


2 - 2cos® - sin 


q 
il 
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Bi + T = 
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Myy = Sq (055 - Cy O53) 
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In these equations, the trigonometric functions are replaced by the appropri- | a 
ate hyperbolic functions for negative values of N. When N vanishes, Sp and eae 


Cm assume the conventional values of 4EI/L and 1/2, respectively. They 


have been tabulated, with slight modifications, in [3] for both tensile and com- 
pressive axial forces. 
E With the joint moments, both in twisting and in bending, thus established in po 
R terms of the joint rotations by means of Eqs. (9), there remain two major 4 
4 considerations to be observed. In the first place, the rotations ( Bij; a ij) of - 
F all the bars joined at panel point i must be compatible with the rigidity of the = 
: joint; in other words, the angles between the individual members must be « 7 
preserved during aay rotation. Secondly, the total of all the bending and : 
twisting moments must be in equilibrium. 7 
To this end, it has been found that the employment of matrix notation leads 
to compactness of both conceptual and computational labor. Consider, for 5 
example, the typical bar connecting joints i and j (see Fig. 2), with its in- nai 1 
clination o jj relative to the positive x-axis defined as shown. All rotations BS: 
and couples are assumed positive as indicated by the arrows, which are se- E 
lected according to the usual right-hand-screw convention. This agrees, it : " 
may be noted, with the previous choice of boundary conditions and, in particu- 
lar, with Eqs. (9). ay | 
If now a moment vector m' anda rotation vector w' are defined by® a 


(11) 


and if, further, a “stiffness matrix” S' anda “carry-over matrix” C' are | 
introduced by means of | 


(12) 


then Eq. (9) may be written 


mj, = -(S'wj, - C'S' w},) (3) 


Let the rotation vector of a joint, if referred to the (x,y) coordinate sys- Ma 

tem, be designated by w, and let similarly the moment vector m refer to 

the same system, or Bt 

; (14) 

j 

that 

; 6. Where possible, the subscripts i and j, which identify the member but do 

? not indicate components, will be dropped if no lack of clarity results from 

this omission. 


T B 
n' = w' = a 
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St 0 Cy 
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(15) 


in which the rotation matrix R and its transpose R* are defined by 


cosa 
sina cosa (16) 


sa 


When Eqs. (15) and (14) are substituted in Eq. (13), the “slope-deflection” 
equations in the (x,y) coordinate system read 


(17) 


s w 5) 
In Eq. (17), S and C are the stiffness and carry-over matrices in the new 
coordinate system, that is, 
(18) 
which are related to the previous (primed) coordinate system by means of the 
transformation equations 

S = R*S'R C = R*C'R (19) 


Hence, by Eqs. (12), (18), and (19), and after some trivial trigonometric 
transformations: 


1/2 (St + Sq) + 1/2 - S,,) cos2a 


- 1/2 (S, - cos2a (20) 


-1/2 (S, - S,) sinaa 


xy yx 


and an equivalent set of equations for the components of C. Since further- 
more R is orthogonal, i.e., RR* =I (the unit matrix), the matrix C'S' itself 
transforms similarly to Eq. (20). Finally, in developing Eq. (17), use was 
made of the fact that a;, = aj; + m and that therefore Rj; = -Rjj; hence, in 
the new coordinate system, the i,j subscripts continue to be interchangeable 
for the stiffness and carry-over matrices. 

In review, the relationships developed above insure that the end rotations 
Bij and 6;; of all bars are compatible with a set of joint rotation vectors 
w,;. The establishment of the equations of equilibrium is now straightforward, 
since Eq. (17) expresses all the bar moment vectors with reference to the 
common (x,y) coordinate system. If n is the number of joints that are free 
to rotate, the following n vector equations (or 2n scalar equations) follow 
therefore immediately: 


(21) 


\ 
4 m' = Rm w' = Rw A 
4 
: 
4 
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where the summation extends over all the bars joined together at i. 
The equations of equilibrium (21) now lead directly to the determinantal 
stability criterion. In fact, define the “total joint stiffness” Kj; by 


(22a) 


and let further 


(if3) (22b) 


Then substitution of Eq. (17) in Eqs. (21) results in the set of linear homo- 


geneous equations 
> Ki = (i=1,2 (23) 


J 


This system has a nontrivial solution if, and only if, the 2n x 2n determinant 
of the coefficients of the joint rotations vanishes. The structure is therefore 
in neutral equilibrium under the smallest load for which 


= 0 (24) 


A Series Criterion of Stability 

Eq. (24) represents a highly involved transcendental equation in terms of 
the external loads; to find its smallest root (that is, the smallest multiplier 
common to all loads) is therefore a task of forbidding numerical complexity. 
In this section, an alternate approach is presented, which constitutes a gen- 
eralization of the “series criterion” of stability, which was first proposed by 
Lundquist [4] inconnection with the stability of frameworks against buckling 
in their own plane. 

Briefly, the suggested method proceeds as follows. Select a joint q (pref- 
erably the one most likely to buckle) and apply at q a unit moment about the 
x-axis. By means of a generalized moment distribution process which will 
be described in detail further below, this moment is distributed among the 
bars joined at q. This introduces unbalanced moments, generally about both 
the x and y axes, at the far ends of these bars. Holding now joint q fixed, 
these moments are all distributed in successive steps until all the other 
joints are balanced. If this process converges, ? the only unbalanced joint will 
be the original joint q, toward which the residual moments ry, and ryx 
about the x- and y-axis, respectively, have been “washed back.” If the same 
procedure is now repeated for an original applied unit moment about the 
y-axis at q, the wash-back moments will be, say, lxy and Tyy: Note that, 
in general, ryy# Tyx: 

If now the “residual matrix” r is defined by 


7. In case of divergence, the structure is known to be unstable; in fact, it is 
unstable even with joint q fixed. 


: 
j 
Ki 843 
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(25) 


then the structure is in general® stable if, and only if, the matrix I- r hasa 
positive trace and determinant; hence 


+ <2 (26) 


(1 ryy)(1 - ryy) - Tyx > 0 


constitutes the desired series criterion of stability. 

Suppose now that, for a given load parameter, the procedure described 
above yields a residual matrix which satisfies conditions (26); the structure 
is then known to be stable for this system of loads, and the process is re- 
peated for a larger load multiplier, until one or the other of the conditions is 
violated. The exact value of the load parameter associated with neutral 
equilibrium can normally be estimated with a fair degree of accuracy through 
interpolation. 

The joint balancing procedure itself represents an extension and generali- 
zation of the familiar moment distribution method; the employment of matrix 
notation makes it possible to describe it in almost identical terms. In fact, 
let the joint i undergo a rotation (vector) ™;, while all other joints remain 
fixed. It then follows, from Eq. (17), that each bar (i,j) exerts a moment 
-S{j;; on this joint, while at the far end of the bars, say at (j), an unbal- 
anced moment -(CS);;W; is introduced. It may be worth repeating, at this 
time, that all the vectors and matrices are expressed in relation to their 
common (x,y) coordinate system. 

If now the rotation wj; is so chosen as to balance the joint, it follows that 
all the new bar moments at i must equal the existing unbalanced moment 
mf}. If, further, the definition of the total joint stiffness K,; is recalled as 
given by Eq. (22a), the required joint rotation must be equal to (K,;)~}m}. 
Hence, whenever an unbalanced moment mf is removed at the joint i, while 
all the other joints remain temporarily fixed, the share taken by each bar at 
i is given by 
= my (27) 
while the new unbalanced moment introduced at the far end j of the bar is 
governed by 


u 
(28) 
The matrix S; (K;;)~2 in Eq. (27) may be interpreted as a “generalized 
distribution factor,” while Cij in Eq. (28), as mentioned before, represents 
a “generalized carry-over factor.” In the problem under discussion, in 
which stability rather than stresses are being investigated, some saving in 


8. A possible, though unlikely, exception occurs when the lowest buckling 


Eqs. (26) constitute a necessary, but not sufficient condition of stability. 


labor can be effected by combining steps (27) and (28)—that is, by transferring 


mode does not involve any joint rotation of q about either axis. In that case, 


= on || 
i 
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unbalanced moments directly from one joint to the other joints. This is 
achieved by letting 


my = -Ty3 my (29a) 


in which the “transfer matrix” Tij is defined by 


Note that this matrix is generally not symmetrical, nor is, in general, Tij 
equal to Tj;. 

The numerical computation of the transfer matrices (two for each mem- 
ber of the framework) and the ensuing moment distribution procedure are 
reasonably simple and can be carried out on the slide rule. After all the 
matrices S and CS have been determined by means of Eqs. (20), the total 
joint stiffnesses K are calculated by adding the corresponding components 
of the stiffness matrices for any one joint. Letting Kyx, Kxy = Kyx, and K 
be the Sonpeneme of K for a given panel point, the inverse (flexibility) 
matrix K~* is given by 


-1 yy 2 

K = l/k - Kyy (30) 


It is presumed here that all total joint stiffness matrices are positive defi- 
nite, since otherwise the structure is known to be unstable, necessitating a 
(probably substantial) reduction in the load factor. Hence, k # 0 and K-1 
exists. 

Let it now be assumed that the two moment distribution procedures have 
been carried out and the residual matrix r has been established as described 
at the beginning of this section. It is evident that if, instead of two unit mo- 
ments, an arbitrary moment m of components m, and my were applied at 
joint q, the residual moment vector would be rm. Since further m = Kw, it 
follows that (I-r)Kw is the actual moment introduced at q. In other words, 
(I-r)K represents the equivalent stiffness of the joint q if all other joints 
are free to rotate. The positive definiteness of this matrix is therefore a 
necessary condition for the stability of the structure. The sufficiency condi- 
tion, which is subject to the exception noted before,8 can be proved by extend- 
ing the argument employed in a similar problem in [5]. 

From this criterion of positive definiteness, conditions (26) can be derived 
directly. Use is made here of the symmetry of (I-r)K, which is a conse- 
quence of Maxwell’s Theorem of Reciprocity and which, incidentally, also 
provides some means of checking the numerical work. It follows that the 
stability of the structure depends on the positiveness of the trace and deter- 
minant of (I-r)K. As for the latter, since K itself has a positive determinant, 
|I-r| must also be positive. In regard to the trace, let the coordinate sys- 
tem be rotated so that K assumes a diagonal form; this is always possible 
on account of its symmetry. By examining the form of (I-r)K in this coordi- 
nate system, and by considering further that the diagonal terms of K are 
both positive, it can readily be verified that (I-r)K has a positive trace if, 
and only if, the trace of (I-r) itself is positive in this, and hence in the 
original, coordinate system. This completes the proof of the inequalities (26). 
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Numerical Example 


The method outlined above is now applied to the truss shown in Fig. 3. All 
joints are connected rigidly and are prevented from lateral sway by means of 
cross bracing (not shown), which, however, does not contribute to their rota- 
tional stiffness. All gusset plates are assumed to be 1/2" thick. 

The purpose of the study is to determine the factor of safety against buck- 
ling—that is, the factor by which the given load must be multiplied to render 
the equilibrium of the structure neutral. In particular, it will be shown that 
although the compression members have been so selected as to make their 
L/r ratios about equal, the danger of buckling out of the plane of the truss is 
considerably greater than that within the plane, at least so far as the struc- 
ture under consideration is concerned. 

As for the stability within the plane of the framework, the investigation 
proceeds along the familiar lines as outlined originally in [4]. The applied 
load W = 50,000 lb was multiplied by several assumed values of the safety 
factor and, for each value, the stiffness and carry-over factors were com- 
puted with the aid of [3]. A unit moment in the plane of the truss was then 
applied at joint B, distributed among the three members joined at B, and, 
with B fixed, the remaining joints were balanced. If the effect of yielding is 
disregarded, it was found that the residual moment at joint B is 0.962 fora 
safety factor of 4.50, and 1.066 for 4.60. Hence, the load corresponding to 
plane instability lies between 225,000 and 230,000 lb. 

The investigation of the lateral stability proceeded along an analogous 
path, as described in the previous section. For an assumed load multiplier, 
the axial bar forces were computed as well as the equivalent forces as given 
in Eqs. (3) and (6). The bending stiffness and carry-over factors S,, and Cy, 
were next determined, as in the case of plane buckling. For the example of 
the double-angle section as used here, a measure of simplification is intro- 
duced by the fac* that the torsional stiffness S; is negligible in comparison 
with Sy. The determination of the S and CS matrices in the (x,y) coordi- 
nate system, as given by Eqs. (20), is simplified accordingly. 

Table I shows a representative portion of the constants involved for an as- 
sumed safety factor of 2.6. The material constants E and G were assumed 
to be equal to 29 x 106 psi and 11.6 x 106 psi, respectively. Since only 
relative stiffnesses are of concern in connection with this study, a number of 
common factors were disregarded in the establishment of S,,. However, it 
is important to retain the correct proportionality in computing S (not needed 
here) from Eqs. (10), particularly since the tables in [3] give values of S,, 
which are only one-fourth of those corresponding to Eqs. (10). 

The total joint stiffness is now determined by adding the corresponding 
components of the S matrices of all the bars connected at a particular joint. 
For example, the stiffness matrix K, of panel point A is obtained as the 
sum Of Sap and Sap; this leads to 


4 
{ 
| 
‘eth 
| 
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0.45 0.45 f 
K, = 
aq 
0.45 1.29 4 
4 
. i which is positive definite and has an inverse [see Eq. (30) | | 
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3.40 -1.29 } 


| -1.29 


The stiffness matrices of all the other joints, and their inverse, are then com- 
puted in the same fashion, 


The final step consists in the establishment of the transfer matrices de- 


fined by Eqs. (29). This is achieved by multiplying the CS matrix of a mem- 
ber, as given in Table I, by the inverse stiffness matrix of the joint from 
which the unbalanced moment is to be transferred. For example, Tap 
represents the effect of the balancing of joint A on the unbalanced moment at 


B; it is given by 
.67 .67|[ 3.40 
.67 1.29 [1.42 0) 
It should be emphasized once again that two such transfer matrices must be 
computed for each bar since, in general, Tap # Tpa. 

Fig. 4 shows all the transfer coefficients ie a safety factor of 2.6; it also 
demonstrates the moment distribution procedure, which in this case turns out 
to be exceptionally simple owing to the symmetry of the structure. A brief 
description at this point, it is felt, may elucidate the simple, but perhaps 
somewhat unfamiliar steps. The introduction and balancing of a unit moment 
about the x-axis at point B produces unbalanced moments about both axes at 
joints A, D, and C, as shown in the first line. When the unbalanced moment 
m, = -0.15 at A is removed, the moments about both axes, which are trans- 
ferred back to B, are each equal to -(1.42)(-0.15) =+0.21; at the same time, 
a moment vector is transferred to joint D whose x-component is zero and 
whose y-component is given -(-0.28)(-0.15) = -0.04. Similarly, the balancing 
of m, = -0.15 at A produces (0,0) at B and (0, +0.04) at D. 


The process is repeated for joint C, with analogous results. The unbal- 
anced moments at D now are -0.29 and 0, respectively, about the x- and 


ASCE 
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STIFFNESS AND CARRY-OVER CONSTANTS FOR A SAFETY FACTOR OF 2.6 a - 
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1G y-axes. The removal of the former introduces an unbalancing effect only at 

fi is joint B, to which an additional moment vector is transferred, whose x-com- 

3 ponent equals -(0.37)(-0.29) =0.11 and whose y-component vanishes. 
7) " With all joints other than B balanced, it is seen that the latter is now sub- 
I jected to a residual unbalanced moment whose components are 0.53 and 0, 
4 i respectively. If the same process is carried out after an initial unit moment 
—o about the y-axis is balanced at B, it can be demonstrated that the residual 
“ moment vector is given by (0, 2.10). Hence, the residual matrix, as defined 
by Eq. (25), is 
© 0.53 
= 
arora ee which is seen to violate the stability criterion (26). In other words, the struc- 
x ture is laterally unstable under a load of 130,000 lb. 
; The same procedure, if carried out for a safety factor of 1.8 (correspond- 
v4 ing to a load of 90,000 Ib), leads to 
Bu 
Es which satisfies both criteria (26). The actual factor of safety of the structure 
ne: is therefore known to lie between 1.8 and 2.6. While its exact value is not im- 
faye Bs portant in connection with the general discussion in this paper, it does seem 
: wg noteworthy that it is very significantly smaller than its counterpart against 
1: plane buckling. 
CONCLUSION 
7 oh As pointed out already in the Introduction, the present theory is based on 
Pe the assumption that, prior to lateral buckling, all members are subjected to 
Ar axial forces through their centroids. What effect, if any, the presence of so- 
f called secondary moments in the plane of the structure has on the condition 
i of lateral stability is difficult to assess without considering a specific ex- 
* $7 ample. It is known, however [6], that structures whose bending moments are 
# statically indeterminate do not ordinarily collapse under loads which are as- 
4 We sociated with incipient lateral instability. The type of framework under con- 
qa sideration represents such an example; hence, it is reasonable to expect the 
7 , actual collapse load to exceed the theoretical buckling load. 
t In any event, the factor of safety which is computed on the basis of the 
iS simplifying assumptica of this paper must represent a lower bound to the 
st factor of safety base. on the collapse load. Indeed, as pointed out in [6], any 
A safety factor based on statically admissible bending moments represents 
sr such a lower bound; obviously, a set of identically vanishing bending moments 
- r is statically consistent in the absence of loads other than those applied at the 
joints. 

[ ¢ Whether the startlingly large discrepancy between the safety factors 
oe against buckling in and out of the plane of the structure can be extended to 
ae trusses of different geometry, or of greater torsional resistance, or 


i 

4 

i 

“we 
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composed of members heavy enough to correspond to plastic buckling, can of 
course not be answered summarily at this point. On the other hand, the re- 
sults of this study seem to indicate that the current design practice of 
equalizing the slenderness ratios of compression members in both directions 
may not necessarily lead to an efficient functioning of the structure. Perhaps 
it may not be unreasonable to suggest that the pertinent features of the pres- 
ent code be subjected to further scrutiny in this respect. 
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WIND INDUCED VIBRATION OF CYLINDRICAL STRUCTURES 


Joseph Penzien,* J. M. ASCE 
(Proc. Paper 1141) 


SYNOPSIS 


This paper presents the results of a wind tunnel investigation to determine 
the dynamic behavior of cylindrical sections under the action of wind at high 
Reynold’s numbers. These results show that wind induced vibration normal 
to the direction of wind is “self excited” rather than “forced”. Two critical 
wind velocities were observed. The lower critical velocity appears to be as- 
sociated with the von Karman vortex street. However, the basic aerodynamic 
phenomenon associated with the higher and more severe critical velocity is 
not presently understood. 


SYMBOLS AND NOTATION 
Definition 


Frequency of vortex shedding 

Wind velocity 

Cylinder diameter 

Strouhal Number 

Unit lateral force resulting from vortex shedding 
Lateral force coefficient 

Time 

Mass density of air--0.00238 slugs 

Maximum dynamic response amplitude 


Static response amplitude due to lateral force at resonance 


Note: Discussion open until June 1, 1957. Paper 1141 is part of the copyrighted 
Journal of the Engineering Mechanics Division of the American Society of Civil 
Engineers, Vol. 83, No. EM 1, January, 1957 


*Asst. Prof. of Civ. Eng., Univ. of California, Berkeley, Calif. 
1141-1 


7 
i} 
EM 1 
| 
4 
Vv 
y Al 
D ‘ 
Cy q 
p 
D ve 
XDR 


1141-2 EM 1 January, 1957 
Symbol Definition 


Reynolds number 

Kinematic viscosity 

Drag coefficient 

Total maximum dynamic reaction of cylinder 
Natural frequency of cylindrical system 


Logarithmic decrement of damping 


INTRODUCTION 


The general problem of designing tall stacks which are dynamically stable 
under the action of wind is becoming increasingly important due to the demand 
by industry for taller stacks. Until recent years, it has been the general 
practice to design stacks conservatively for a static load equal to the wind 
drag at a wind velocity considera ily above the maximum velocity expected. 
However numerous cases(1,2,3,4,9)* of serious wind-induced vibration have 
been experienced and in some cases failures have occurred, thus pointing out 
the inadequacy of designing stacks on this basis. 

Numerous papers have been published on this subject recently which point 
out the tendency of stacks to vibrate in a direction normal to the wind. These 
papers usually associate such behavior with the well known von Karmon vor- 
tex street.(6) This phenomenon is a periodic shedding of vortices alternately 
from the cylinder at a frequency (f,) as given by Eq. 1. 


f,+5(%) (1) 


This periodic shedding of vortices produces a change in the particle ve- 
locities on opposite sides of the cylinder and by Bernoulli’s theorem similar 
changes in pressures result. Thus, it is generally assumed an alternating 
unit force (F) normal to the direction of the wind is produced as described 
by Eq. 2. 


2 
es sin (2) 


This periodic forcing function gives rise to forced vibration, the amplitude of 
which varies with wind velocity as shown in Fig. 1 and the frequency of which 
is equal to that of the forcing function itself. The maximum response with 
this type of excitation is reached when the vortex shedding frequency is equal 
to the natural frequency of the stack. 

Much of the discussion in the literature on this problem has been con- 
cerned with the numerical values which should be assigned to the Strouhal 
Number S and the Lateral Force Coefficient C, . 


*Numbers in parentheses refer to the bibliography. 
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It has been concluded by some writers that above the critical range of 


Reynold’s numbers (R -“?), the shedding of vortices become aperiodic. In 


view of this one would conclude that a truly resonant condition would never be 
reached and excessive vibrations would not build up. Experience however 
shows that excessive vibrations can occur at Reynold’s numbers considerably 
above the critical range. 

Recently a paper(3 by Ozker and Smith has presented a new concept of the 
dynamic response of tall stacks to wind action. By measuring the dynamic 
response of a full scale stack to a wide range of wind velocities they con- 
cluded that wind induced vibrations of a stack are “self excited” rather than 
“forced” as discussed above. With this type of oscillation the frequency of 
vortex shedding is controlled by the response of the stack itself. That is, the 
stack vibrates at its own natural frequency at all times and is independent of 
wind velocity. Similar response measurements were made by Amori(4) in 
1918 which support the above concept. 

Presented in this paper are the results of a wind tunnel investigation of 
the dynamic behavior of cylindrical sections under the action of wind. In 
general these data support the concept of self excitation set forth by Ozker 
and Smith. There are however differences in the response—wind velocity 
relationship. 


Experimental Techniques 


Three aluminum cylinders 44 inches long and having diameters of 3, 6, 
and 9 inches were tested in a 3 foot open throat wind tunnel. Each cylinder 
was supported by a cantileverbeam at each end as shown in Fig. 2. These 
beams were placed parallel to the direction of the wind and outside the main 
air stream. In this manner the rigidity of the test specimen was kept high in 
the direction of the wind while its rigidity normal to the wind was kept rela- 
tively low. By changing the width and length of these supporting beams a 
wide range of frequencies could be obtained for each cylinder. 

The dynamic response was measured by SR-4 gages mounted on each side 
of the supporting beams. These gages were wired into a bridge in such a way 
that only a transverse displacement of the cylinder as a whole would be re- 
corded, i.e., a rotation of the cylinder about an axis parallel to and midway 
between the supporting beams would not change the bridge balance. For- 
tunately, however, because of symmetry of structure and loading, rotations 
of this type were negligible. The strain gage bridge was connected to a 
Sanborn Recorder which produced a continuous record of dynamic response 
for a given system. 

A total of twenty tests were conducted as shown in Table 1. Each test 
represents a given cylinder with a given support rigidity. A continuous 
record of dynamic response was recorded for each test during which the wind 
velocity was increased in increments of 5 ft/sec from zero to 115 ft/sec. 

Vertical fins were attached to the 6 inch cylinder for three tests. These 
fins projected 5/8 of an inch from the cylinder and were spaced every 30 de- 
grees around its perifery as shown in Fig. 3. These tests were conducted to 
check the effect of fluting on dynamic response. 
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TABLE 1 - TEST PROGRAM 


Test Cylinder Natural Test Cylinder Natural 
No. Diameter Frequency No. Diameter Frequency 
In. Cycles/ Sec. In. Cycles/ Sec. 


oO O 


10 


*5/8 inch fins attached every 30 degrees around perifery of cylinder. 


DISCUSSION OF TEST RESULTS 


Type of Response 


It was clearly shown in this series of tests that the response frequency of 
a cylindrical system under the action of wind remains constant regardless of 
wind velocity and is equal to the natural frequency of the system. This indi- 
cates that the vibration of tall stacks are “self excited” rather than *forced” 
as generally assumed. 

A typical record of cylinder response normal to the direction of wind is 
shown in Fig. 4. This figure shows a portion of the continuous record of 
response for the 3 inch cylinder (f, = 17 c/s) subjected to a steady wind 
velocity of 65 ft/sec. 

In addition to the fact that each cylinder responded at its natural frequency, 
another important fact is immediately apparent from these records. That is, 
the vibration amplitude is generally not constant during a condition of steady 
wind as one might expect. Instead the oscillation increases rapidly and then 
decreases in a similar manner. While this behavior is somewhat random, it 
does indicate that energy is suddenly fed into the system from the airstream. 
Some of this energy is dissipated through structural and aerodynamic damp- 
ing, while the remainder is suddenly fed back into the airstream. 

The basic aerodynamic phenomenon causing the above alternate conditions 
of stability and instability is unknown at present. One might set forth several 
hypotheses as possibilities; however there is no experimental or theoretical 
evidence to support them. One possibility is that there is sufficient periodic 
shedding of vortices on alternate sides of the cylinder to build up the vibra- 
tion. After a number of cycles the shedding becomes out of phase with 
response, thus reducing the vibration. Another possibility is that at these 
high Reynold’s numbers an alternating force is produced due to flow changes 
resulting from an alternate shifting of the separation point on the boundary 
layer. Shedding of vortices in this case would be considered aperiodic. This 
shifting of separation point can be alternately in and out of phase with re- 
sponse. 
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It should be noted that the type of response reported above is very similar 
to the full scale stack response as reported by Ozker and Smith. This indi- 
cates that the aerodynamic phenomenon existing for these tests is similar to 
that for full scale stacks. One difference to be questioned however in cor- 
relating these data with full scale is that while these data are in the critical 
range of Reynold’s numbers (105 7 106), full scale stacks usually vibrate at 
Reynold’s numbers above the critical range (> 10%). 


Strouhal Number‘”) 


The relationship between Strouhal number and Reynold’s number has been 
measured experimentally by Relf and Simmons'®) in the critical range for 
rigid cylinders. This relationship is shown in Fig. 5 along with the relation- 
ship between drag coefficient and Reynold’s number. The change from 
“laminar” to “turbulent” flow in the critical range results in a sudden reduc- 
tion of drag coefficient and an increase in Strouhal number. 

Since the response as reported above for non-rigid cylinders is “self ex- 
cited,” the Strouhal number might appear to be of little significance. How- 
ever, when considering the rapid increase and decrease of response as de- 
scribed above, one might expect the periodic force to remain in resonance 
through more and more cycles as the wind velocity gets nearer and nearer to 
the resonant wind velocity as given by Equation 1. This condition would result 
in a build up of response at this wind velocity. The value of Strouhal number 
used in calculating this resonant wind velocity would be that obtained experi- 
mentally for rigid cylinders. 

As a check on Strouhal frequency, a test was conducted in this investiga- 
tion for which a cylinder was made sufficiently rigid so that it responded not 
at its natural frequency but at a frequency assumed to be the Strouhal fre- 
quency. The 9 inch cylinder which was used for this test had a natural fre- 
quency of 43 cycles/second. A portion of the response record for a wind 
velocity of 35 ft/sec is shown in Fig. 6. This response record indicates a 
somewhat steady condition. Therefore, one can get a measure of the pre- 
dominant frequency of the forcing function. 

These frequencies as obtained from this test are plotted in Fig. 5. While 
these results deviate somewhat from those of Relf and Simmons, they do in- 
dicate a sudden increase of Strouhal number in the critical range. 

Also plotted in Fig. 5 are the values of Strouhal number as calculated 
using Equation 1 and the lowest wind velocity at which a resonant condition 
was reached as shown by Point No. 1 on the response curves in Figs. 7, 8, 
and 9. The correlation between the values of Strouhal number obtained in 
this manner and those obtained from the rigid cylinder test described above 
is reasonably good. This indicates that a first resonant point will be reached 
at the resonant wind velocity as associated with the von Karmon vortex street. 


Amplitude of Response 


The experimental relationship between cylinder response and steady wind 
velocity is given in Figs. 7, 8, and 9 for the 3, 6, and 9 inch smooth cylinders, 
respectively. Each curve represents a different degree of stiffness for the 
system as indicated by the given natural frequency f,. The amplitude of 
response is shown by plotting the maximum total dynamic reaction force (Fp) 
normal to the direction of the wind for each increment of wind velocity tested. 
This force is equivalent to that static force which would be required to 
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produce the same stress as the maximum dynamic stress. 
Although, it is difficult to find correlation among these data, the following 
important facts are indicated by these results. 


(1) The wind velocity at which the first peak response occurs (Pt. No. 1, 
Figs. 7, 8, and 9) agrees reasonably well with the resonant wind velocity as 
calculated using Equation 1 provided allowances are made for the change of 
Strouhal number with Reynold’s number as previously discussed. 

(2) A more severe resonant condition may be reached at a higher wind 
velocity than in (1) above. This condition is most clearly shown in Fig. 8 at 
a wind velocity of approximately 60 ft/sec. The vibration became so severe 
in this velocity region that it had to be stopped to prevent yielding of the sup- 
porting beams. 

(3) Maximum response in general does not continue to increase with wind 
velocity, but does reach a maximum at certain critical velocities. 

(4) There appears to be a certain critical stiffness or rigidity for which 
vibration due to wind is most severe. 

(5) Stresses produced by vibration of a cylinder normal to the direction of 
the wind can be much larger than those produced by the steady drag force 
provided the right combination of wind velocity, structural stiffness, and 
structural geometry is present. Low structural damping, as is usually the 
case in tall steel stacks, is also necessary for this condition to exist. 


Effect of Surface Roughness on Dynamic Response 


Response curves for the fluted 6 inch cylinder (Fig. 3) are shown in Fig. 
10. A comparison of these data with those for the smooth 6 inch cylinder 
shows an appreciable reduction of vibration. This fact indicates that while 
roughness or projections from a smooth cylinder increases the steady drag 
load it does have a stabilizing effect on vibration normal to the wind. 


GENERAL CONCLUSIONS AND RECOMMENDATIONS 


From the results of this investigation, the following general conclusions 
may be deduced: 


(1) Wind induced vibrations of cylindrical sections normal to the wind are 
“self excited” rather than “forced”. 

(2) These wind induced vibrations build up rapidly at certain critical wind 
velocities. 

(3) A lower critical wind velocity exists which corresponds to that asso- 
ciated with the von Karman vortex street. The oscillation which builds 
up at this velocity is usually not too severe. 

(4) A high critical wind velocity can exist which may produce very severe 
oscillations. The basic aerodynamic phenomenon causing this behavior 
is not presently understood. 


In view of the present state of knowledge regarding wind induced vibration 
of tall stacks, it is very desireable that further research be conducted to 
clarify the basic aerodynamic phenomenon associated with this behavior. 

Since it is very difficult to conduct wind tunnel studies above the critical 
range of Reynold’s number, it will be necessary that full scale tests be made. 
These tests should be so designed to provide a clear understanding of the 
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aerodynamic forces involved. Only after such studies are successfully com- 
pleted will the designer have the necessary information to arrive at economi- 
cal and safe designs of tall stacks. 
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Fig.f- Dynamic Response Resulting trom Periodic 
Forcing Function as Described by Equation 2. 
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Fig. 3 - Test Section of Wind Tunnel with 6 inch fluted cylinder in place. 
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Fig. 6 - Dynamic Response of Rigid Cylinder. 
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Discussion of 
“LATERAL BUCKLING OF ROLLED STEEL BEAMS” 


by R. A. Hechtman, J. S. Hattrup, E. F. Styer, and J. L. Tiedemann 
(Proc. Paper 797) 


CORRECTIONS.—In Figs. 9 and 10 curves that are presently labeled 
“From Eq. 7” should be labeled “From Eq. 2.” Also, Eq. 2 should be 
changed so that the term L® is replaced by L2. 
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Discussion of 
“SOME OBSERVATIONS ON OPEN CHANNEL FLOW AT 
SMALL REYNOLDS NUMBERS” 


by Lorenz G. Straub, Edward Silberman, and Herbert C. Nelson 
(Proc. Paper 1031) 


WALLACE M. LANSFORD,* M. ASCE, and JAMES M. ROBERTSON, ** 

M. ASCE.—This very concise presentation of an extensive study of laminar 
and turbulent flow in open channels is most welcome. Since the writers’ in- 
terest in open channel flow has primarily been in connection with the turbu- 
lent flow regime, their discussion of the laminar flow results presented will 
be limited. Comments on the turbulent flow aspects will be more extensive. 
With regard to the theoretical relations presented for viscous laminar flow in 
elliptical, rectangular and 90 degree triangular open channels, it should be 
noted that these are contained in the theoretical results presented by 
Boussinesq in 1868(1) for viscous flow in rectangular and elliptical conduits. 
In applying these to the appropriate open-channel cases, one assumes no ef- 
fect of the free surface. In turbulent flow, this assumption is not quite cor- 
rect and the maximum velocity usually lies below rather than on the free 
surface. 

Did the authors’ measurements of the velocity distribution in the case of 
laminar flow include a determination of the location of the filament of maxi- 
mum velocity? It would be most interesting to know if this occurs on, or be- 
low, the free surface. For turbulent flows, the occurance of the maximum 
filament below the surface has been variously attributed to air resistance 
surface tension, secondary currents or some effect of turbulence. Observa- 
tions of the occurance or non-occurance of this phenomenon in laminar flow 
would at least represent the case of no turbulence. They would tell us about 
the importance of turbulence and could well shed some light on this intruigu- 
ing but elusive question. 

Studies of open-channel flow, although not as extensive as those reported 
by the authors, have been under way in the Fluid Mechanics and Hydraulics 
Laboratory of the University of Dlinois for a number of years.(2,3) These 
have been mostly concerned with the turbulent flow regime at mild slopes. 
The following presentation of the Illinois’ data is included as it directly sup- 
plements the authors’ work. 

The mean flow relations in a smooth 90 degree triangular flume at 
Reynolds numbers up to 45,000 are presented by the authors in Figs. 2 and 5. 
The results obtained by W. M. Owen(2) on a similar smooth wooden channel 
at Reynolds numbers from 20,000 to 275,000 are presented in Fig. 1 in com- 
parison with the Minnesota data and smooth pipe relations. Baring any sys- 
tematic errors, it would appear that at the higher Reynolds numbers the 
friction for triangular channels falls below the Blasius relation for smooth 
pipes, as given by Eq. (9a). 


~ * Prof. of Theoretical and Applied Mechanics, Univ. of Dlinois, Urbana, Il. 


** Prof. of Theoretical and Applied Mechanics, Univ. of Ilinois, Urbana, Il. 
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The open channel flow studies at Illinois have also been concerned with the 
nature of the turbulent flow in the triangular channel roughened with sand 
grains. This type of roughness was chosen because it more nearly simulates 
the roughness type about which we know so much in pipe flow, i.e., 
Nikuradse’s sand grain roughness.(4) Figure 1 also includes a comparison of 
the pipe and open channel variation in friction factor with Reynolds number 
and relative roughness (k/2R). In contradiction to Owen’s(2) conclusion, from 
some of the same data, these results show a transition region between the 
smooth and fully rough conditions for open channel flow. Differences from 
the pipe results are apparent. Thus, the present sand grains appear to be 
effectively rougher than Nikuradse’s (by about a factor of 1.4) and fully rough 
flow seems to occur at a lower Reynolds number. The first of these observa- 
tions is not particularly unexpected since the sand grains used in the flume 
were not as closely packed as Nikuradse’s were. 

A better understanding of the nature of the flow is obtained through studies 
of the velocity distribution. Thus, in Fig. 6, the authors present velocity 
traverses obtained for turbulent flow in the smooth 90 degree flume. Similar 
measurements in the roughened triangular flume have been taken at Illinois 
by D. B. Y. Chen.(5) These measurements were made in the same 55 foot 
flume used by Owen.(2) The sand grain roughness was not identical, but the 
same sized grains were attached in the same fashion. That the roughness 
was effectively the same is evident from the general agreement of friction 
factors in Fig. 1. More baffling was added at the entrance of the flume to im- 
prove flow symmetry. The velocity measurements were taken at a station 45 
feet from the head bay, with a pitot-static tube (diameter 0.106 inches). For 
a flow depth of 0.518 feet, velocities were measured at 32 locations; while for 
a depth of 0.788 feet, 48 locations were used. The reliability of the velocity 
measurements was verified through comparison of integrated discharges with 
those obtained by weight. The integrated discharges were quite consistently 
smaller by about 2.75 percent. 

The results obtained for the smaller depth and roughness (2R/, = 294) are 
shown as contour plots in Fig. 2. The symbols on the upper plot represent 
the locations at which velocities were actually obtained. The differences be- 
tween the several rlots are most apparent when one considers the maximum 
velocity in magnitude and location. It is seen from this figure that as the 
Reynolds number increases, the ratio of the maximum velocity to the average 
velocity increases and this maximum is found nearer the free surface. All 
thirty traverses taken were analyzed by noting the ratio of maximum to aver- 
age velocity Wm/, and the ratio of the depth of the maximum velocity location 
to depth of flow €. The results are summarized in Table I and, in general, 
agree with trend shown in Fig. 2. However, the measurements were mostly 
taken in the region of rough conduit where the friction factor presumably de- 
pends only upon relative roughness and not Reynolds number. Variation of 
the velocity ratio with the relative roughness, as found by Nikuradse in rough 
pipes, is not evident from these data due either to lack of precision or effect 
of the free surface. 

Five of the velocity traverses taken by Chen were analyzed to determine 
velocity profiles in the rough flume in a manner similar to those presented 
by the authors in Fig. 6 for the smooth condition. Profiles have been ob- 
tained on the center line and on normals to the side walls passing through the 
point of maximum velocity. In this latter case, values of the velocity were 
determined from velocity contours sketched with the aid of the measured data. 
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The resultant velocity profiles are presented in Fig. 3. Comparison of the 
profiles on the center line and the normal lines indicates similar effects to 
those shown in the authors’ Fig. 6. The mean lines passing through the data : 7 
have the following equations: r 

for the normal: 


= 8.3 + 5.75 log (*/,) 


for the center line: 


= -4.85 + 10.5 log 


This form of equation, rather than an exponential form as used by the authors, 
is chosen to permit comparison with Nikuradse’s pipe work. In making this 
comparison, the first equation above is found to be very close to Nikuradse’s 
equation for the velocity profile in fully rough pipe flow. Such similarity is 
coincidental since it has already been noted that the roughness heights of the Pes sl py 
sand grains used correspond to higher values of Nikuradse(6) sand grain 
roughness. Another factor influencing these results is the variation in wall 
shear stress around the periphery of the flume. Thus, for the center line 
traverses, the effective shear must be much less than the mean value used to 
compute w,. 

Some idea as to the variation in the wall shear stress around the wetted 
perimeter of the flume has been obtained from the velocity measurements 
nearest the side walls. For each such measurement, the Nikuradse rough 
pipe velocity profile equation,* corrected for the different effective roughness 
height of the sand grains, was utilized to calculate the value of the local shear 
stress velocity w,. The square of the ratio of w, to the mean shear stress 
velocity w, yields the ratio of wall shear stress to the average or T,,/7RS. 
The indicated variation in wall shear stress up the side of the flume is pre- 
sented in Fig. 4. Although there is considerable scatter and traverse A ap- at 
pears out of line, a definite trend is apparent for the wall shear stress varia- A ie 
tion. 

For rectangular rough channels, Keylegan(6) has analyzed Bazin’s data in ' 
a similar manner to obtain the peripheral variation in wall shear. He notes ¢ 
that the shear has a maximum value at the midpoint of the bottom and minima , 
at the corners. With this information, he then verifies the agreement between 
Bazin’s velocity traverses and the Nikuradse velocity equation. It is the 
writers’ opinion that more such studies of shear and velocity distributions in 
turbulent open channel flow will advance our basic understanding of the 
phenomena involved. 

It is hoped that this presentation and discussion will aid in giving a clearer 
insight into some phases of open channel flow, thus supplementing the excel- 
lent coverage of the author’s paper. - 


* From a detailed analysis of Nikuradse’s velocity measurements near the 
wall, this equation has been found to apply for values of w,k/y greater than 
4.4. For the traverses analyzed, W,k/v was not less than 8.4. This im- 
plies an error when w/w, is less than 0.52 (ry,/»RS less than 0.23). \ 
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TABLE I 
Summary of Maximum Velocity Data 
From Triangular Flume With Rough Sides 

d k Ww 2R kW, Wm 

0.788 0.00125 0,0005 1.293 444 144000 84 1.164 0.238 
0.788 0.00125 0,001 1785 444 202000 118 1,240 0, 195 
0.788 0.00125 0.002 2,56 444 290000 17.0 1.268 0.031 
0.518 0.00125 0,0005 0,989 294 74360 6.9 1175 0,212 
0.518 0.00125 0,001 1398 294 104000 9,8 1.204 0, 169 
0.518 0.00125 0,002 1988 294 146000 13.6 1,232 0, 135 
0.788 0.0025 0.0005 1.235 206 135600 16.5 1.184 0,290 
0.788 0.0025 0, OOL 1747 206 190000 30,8 1,264 0, 197 
0.788 0.0025 0.002 2,431 206 287000 35.5 L211 0.042 
0.518 0.0025 0.0005 0.908 134 65200 13.5 1231 0,254 
0.518 0.0025 0. 001 1.275 134 94000 19,1 1.181 0,212 
0.518 0.0025 0, 002 1, 863 134 134000 26.6 1,240 0,171 

W = average velocity, fps 
*. = maximum velocity, fps 

€ = depth to maximum velocity filament divided by depth of 


flow (d) 
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W=0.993 
S =0.0005 
Re = 74 700 


W=1.388 
S =0.00! 
Re =103 200 


W=1.996 
#=0.002 


Re=146 600 


d=O.5I8 FT k = 0.00125 FT 


FIG.2 VELOCITY CONTOUR PLOTS 
FOR THREE REYNOLDS NUMBERS 
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RAVERSE | SYMBOL 
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| FIG.3 VELOCITY PROFILES FOR ROUGH 
z: TURBULENT FLOW IN TRIANGULAR FLUME 
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TRAVERSE| A | B | 
0.518' 
0.001 
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FIG.4 SHEAR STRESS VARIATION UP THE 
SIDE OF 90-DEGREE FLUME 
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C. 4, POsEY,! M. ASCE.—The writer questions part of one of the conclu- 
sions, in which it is stated that “The entire effect of shape in laminar flow 
may be computed theoretically.” If it is meant to imply that the theory based 
on the Navier-Stokes equations is applicable to the extent that it will yield an 
accurate picture of the velocity distribution, there is good evidence to the 
contrary. 

The author’s data shown in Fig. 3 show a considerable deviation from the 
theoretical which is evidently systematic, and not haphazard. Data obtained 
by Ellis B. Pickett show deviations of comparable magnitude. 2 Theoretical 
studies based on differential equations that describe the behavior of actual 
fluids more precisely than do the Navier-Stokes equations yield significantly 
different results.3,4,5 

It is reasonable to expect that velocity distributions which differ appreci- 
ably from geometrical similarity might still have very nearly the same over- 
all rate of friction loss. This could explain the agreement shown in the 
author’s Fig. 2. 
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Re Discussion of 
oo) “LATERAL BRACING FORCES ON BEAMS AND COLUMNS” 


by William Zuk 
(Proc. Paper 1032) 


my ROBERT V. HAUER.!—The theory of elastic stability dues not establish a 
pe requirement of strength for a member bracing a column or a beam which is 
. free from imperfections and loaded concentrically. In order to provide the 
designer with any information of this sort, the author enforces a solution by 
assuming an initial crookedness of the member to be braced. 

fo It seems doubtful to the writer whether much is gained by such a proce- 
Ld dure, as its net result is only to replace the usual thumb rule—namely, to 

5 assume a bracing force equal to a certain percentage of the force in the main 
member—with another thumb rule assumption about the magnitude of the 
initial crookedness. Such an accidental deviation from the straight line is 
only one of the many differences between the actual physical conditions and 
the idealized assumptions upon which any structural design is necessarily 
based, and, therefore, should be taken care of by a properly chosen safety 
factor. 

Furthermore, the artificial approach which the author chooses for the 
problem in question is not necessary at all. For, the theory of elastic stabili- 
ty does furnish the designer with some requirement regarding a bracing 
member. This requirement is not one of a certain strength, but one of a 
minimum stiffness. 

For instance: If a column with a total length L is braced at mid-height, an 
investigation of the elastic line shows that the unsupported length can be as- 
sumed to be 1/2 L only under the condition that the spring constant of the 
bracing member is at least equal to 4 P,/L, where Pz is the critical load 
under which the braced column loses its stability. If the spring constant is 
less, the critical load, i.e., the capacity of the column, becomes smaller (but 
it does not become greater, if the spring constant increases). Similar re- 
sults are obtained in all other cases of bracing. 

Generally, the stiffness requirement for the bracing member is not hard 
to fulfill, as the following numerical example shows: A 10 WF 33 column with 
a total length of 28 feet is divided into two equal lengths of 14 feet by a brac- 
ing member. Assuming the steel having a modulus of elasticity of 30,000 kips 
per sq. in., a yield point of 40 kips per sq. in. and a proportional limit of 27 
kips per sq. in., there exists in the inelastic range the following approximate 
relation between the slenderness ratio 1/r and the critical stress fo:2 


f. (in kips per sq. in.) = 40 - 544 (l/r) 


This gives for the example in consideration f, = 31 kips per sq. in. or P = 
31 x 9.71 = 301 kips. 


1. Structural Engr., Albert Kahn, Inc., Detroit, Mich. 
2. Friedrich Bleich: Buckling Strength of Metal Structures, McGraw-Hill 
Book Company, Inc., New York, Page 53, Equ. 63. 
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In order that this computation be correct, the spring constant of the brac- 
ing member must be at least 


ec 4x 301 
1 L 28 


= 43 kips per ft., 


or 3.58 kips per in. In other words, under an axial load on the brace of only 
3.58 kips the bracing member is allowed to give as much as a full inch, which 
shows that even a rather flexible lateral support is sufficient to brace the 
column properly. 

In more complicated cases, as in a half through bridge, where the top 
chord is braced by open frames, the influence of the flexibility of the lateral 
supports upon the capacity of the top chord must be carefully investigated, 
regardless of any more or less fictitious loads which are assumed to act 
horizontally on these supports. For, a requirement of stiffness can never be 
replaced by a requirement of strength. 


CHESLEY J. Posey, M. ASCE.—The writer wonders whether it is ad- 
visable to attempt to specify the force that lateral bracing must supply. In 
Fig. A, consider column AB which has frictionless hinges at each end and 
also a hinge at its midpoint, C. Points A and C are restrained from moving 
in the horizontal direction by braces DA and EC. Assuming for the time be- 
ing that the dead load of member EC can be neglected, how should it be de- 
signed? It is seen immediately that the stress in EC given by the ordinary 
analysis of statics is zero, but that the member cannot be dispensed with, for 
then the column will buckle at the hinge C, and collapse. To say that the 
member must be strong enough to carry its own weight and conform to the 
specification requirements for minimum slenderness ratio is begging the 
question. 

Of course fabrication and erection difficulties would not permit points A, 
C, and B to lie exactly on a straight line. Assume that point C is out of line 
by a small distance §. A triangle-of-force solution applied to the free-body 
diagram for joint C shows that member C must now supply a force 


2PS5 (a) 


s 


F- 


This relationship, graphed in Fig. B, is approximate in that it will only hold 
as long as § is small compared with f. However, lack of linearity should not 
be noticeable until 6/f is 0.05 or more, a value much greater than is likely 
to occur. When lack of linearity does appear, however, it will require a 
larger value of F than that given by Eq. (a). 

Fig. B shows the total force that member EC must provide, depending upon 
the amount of slight misalignment of the joints. Without making an assump- 
tion as to how large the misalignment might be, it is impossible to fix the 
magnitude of the force that member EC must be designed for. 

If member EC is subjected to a force, however, it will deflect. The rela- 
tionship between F"', the total stress in EC, and the deflection due to that 
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F' is the force that member EC can supply when joint C is moveda distance 
6 from its unstressed position. If now we design member EC so that 


it will be just able to supply the force necessary to hold the joint C in 
equilibrium for any small value of 6. But 6 can not be permitted to become 
large, for the force F required to prevent buckling will eventually become 
larger than the value given by the linear equation, while the force F" that 
member EC can provide will become smaller than the linear value as soon 
as the proportional limit is exceeded. Evidently this design is not safe, and 
some margin of safety must be provided. The design requirement for mem- 
ber EC is seen not to be a strength requirement, but rather a “stiffness” re- 
quirement, since the factor which must be safely greater than 2P/s is the 
ratio of the axial force to axial deflection in member EC. This axial stiffness 
would be slightly different in tension and compression. It is the smaller 
value that must be adequate. 

For the example of the author’s Fig. 10, the writer agrees that the mini- 
mum value of the stiffness coefficient of the elastic support is P/f, and that 
a limiting value 50% greater is practical. This is a simple and basic require- 
ment which could be made the basis of specifications, and the writer sees no 
reason to shift to a strength requirement which has no validity unless defi- 
nitely related to the specification governing the permissible initial eccen- 
tricity. I the fabrication alignment does not meet specifications, the margin 
of safety of the strength requirement is encroached upon. 

It may be argued that a lateral bracing specification based solely upon 
stiffness is inadequate because failure might occur before the force reached 
the value demanded by the initial eccentricity. That such an extreme brittle- 
ness could occur with ordinary structural materials seems extremely unlike- 
ly, but if it is thought to be possible, safety could be assured by requiring that 
the necessary stiffness be maintained for loads up to say 200% of the value 
demanded by the initial eccentricity. 


C. M. SEGEDIN.*—There are one or two points that can be made in con- 
nection with the author’s Case 1, in which he considers a column with initial 
crookedness given by y, = asin(7x/L). He gives in equation (1.1) the lateral 
bracing force, F, to be applied at the midpoint of the column if there is to be 
no further deflection there when an axial load P is applied. If this expression 
is simplified and rearranged somewhat, it is found that F can be given in the 
form 


f(u) 


f(u) 
where ujs 

3(1-4n2/n?)( tanu-u ) 
and ustpl, = PAI, 


a. Senior Lecturer in Math., School of Eng., University College, Auckland, 
N. y 
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f(u) of course depends, non-dimensionally, on P, EI and L, but evaluation 
yields the result that it is a steadily increasing function of u - at least until 
tanu = u - and it is found that 


£(0) = 1, = 24/96 = 1.015, f(x) = 27/9 = 1.097. 


These values of u correspond to the cases (a) no axial load, (b) axial load 
equal to the first Euler load, P = EIn2/L2, and (c) axial load equal to the 
second Euler load, P = 4EIn2/12, at which the column in this case buckles. 
Thus for the range of practicable values of P, f(u) is remarkably constant and 
we have the result that F is proportional to a, as pointed out by the author, 
and approximately proportional to P and 1/L. 

Increased stability may be achieved if the slope of the column at the mid- 
point is kept zero. In this case buckling occurs when 


tanu = u, namely, u = 1.430K or P = 8.18xEI/L*, 


and it is easily shown that no torque need be applied at the midpoint to pre- 
vent the tendency of the column to turn there. 

This suggests a related problem for which the buckling load is again 
P = 8.1872EI/L2 in which the column, simply supported at each end, is as- 
sumed to have an initial crookedness given now by y, = asin(27x/L). If fur- 
ther deflection and slope are both kept zero it can be shown that the bracing 
torque, M, applied at the midpoint must be 


M = (6Pa/x)F(u) 


where F(u) = —_u*tem 
3(1-u2/x2)(tanu-u ) 


F(u) is not so conveniently constant for different values of u but here M is 
relatively insensitive to L at least for small axial loads. Throughout, for the 
initial crookedness now assumed the lateral bracing force is zero. 
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Discussion of 
“MINIMUM-WEIGHT DESIGN OF A PORTAL FRAME” 


by William Prager 
(Proc. Paper 1073) 


ROBERT L. KETTER,* J.M. ASCE.--In discussing the problem of the 
minimum-weight design of portal frames, it should be pointed out that seldom 
if ever will a “real” structure be subjected to only one condition of loading. 

It is therefore necessary to determine (or check) for each of the possible 
loading configurations the required sized of the component parts of the frame. 
The problem becomes even more involved for the case of multiple span 
structures. 

As a second and final point of this discussion, it is desirable to show how 
well the @ = 2/3 used by the author checks currently available structural 
sections. In Figure 1 is plotted to a log-log scale the weight per foot of 
members versus the plastic modulus (which is proportional to the full plastic 
moment, M,) for each of the rolled WF shapes. As indicated, within a given 
range of nominal depth of section, weight increases approximately as the 
plastic modulus raised to the 0.91 power. Considering the lower extremities 
of these sections available (that is, the most economical sections), if the 
weight is less than about 40 lbs. the weight increases approximately as the 
0.55 power of the plastic modulus. For those cases where the weight is 


greater than 40 lbs., W = k(Mp)9-63 seems to be a better approximation. The 
a= 2/3 value used by the author is a reasonable value for the overall range of 
available wide flange shapes. 


* Research Asst. Prof., Fritz Eng. Lab., Lehigh Univ., Bethlehem, Pa. 
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Journal of the 
ENGINEERING MECHANICS DIVISION 


Proceedings of the American Society of Civil Engineers 


A PRESSURE LINE CONCEPT FOR INELASTIC BENDING 


Frank Baron, M. ASCE 
(Proc. Paper 1157) 


SYNOPSIS 


The concept of the pressure line, useful for arches and rigid frames, is 
developed for quickly estimating the effects of plasticity on the behavior of 
structural elements subjected to axial and flexural loads. An initial estimate, 
obtained by means of an elementary theory of mechanics, is adjusted to fit 
the conditions of plasticity. 


A procedure is presented for quickly estimating the effects of non-linear 
stress-strain relationships on the behavior of structural elements subjected 
to axial and flexural loads. The procedure consists of obtaining an initial 
estimate by means of an elementary theory of mechanics and adjusting the 
estimate to fit the conditions of a theory of plasticity. In this way, the re- 
sults of both theories are available for comparisons and for further studies 
as in probléms of ultimate load analyses and inelastic buckling. 

Figures 1 and 2 show a cross-section of a differential element subjected 
to an axial load, P, and a moment, M, about a centroidal axis of the section. 
For convenience of this discussion, the cross-section is assumed to have an 
axis of symmetry as shown. The effects of a shear, Poisson’s ratio, residual 
strains, and local instability of projecting elements of the cross-section are 
herein neglected. 

The elementary theory of mechanics is summarized in Figure 1 wherein 
the stress in a fiber is defined as 


8, = 8, + 8, (1) 


and where sg = P/A and Sp = Mx;,/I. The distribution of stress across the 
section is linear and satisfies the requirements of statics; namely, 


Note: Discussion open until June 1, 1957, Paper 1157 is part of the copyrighted Journal 
of the Engineering Mechanics Division of the American Society of Civil Engineers, 
Vol, 83, No, EM 1, January, 1957. 


1. Prof. of Civ. Eng. Univ. of California, Berkeley, Calif. 
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= P (2) 
= M (3) 


Among the assumptions of the elementary theory are a linear distribution of 
strains across the section and a linear stress-strain relationship for the ma- 
terial. The angle-change between adjacent cross-sections is stated in terms 
of strains, stresses, and moment as follows: 


where A€ designates the difference in the strains of the extreme fibers, Asp 
the difference in the stresses of the extreme fibers, and E, the modulus of 
elasticity. 

In the theory of plasticity discussed herein, the requirements of statics 
and of geometry are the same as those in the elementary theory of mechanics. 
In each theory the distribution of stresses across the section must satisfy the 
requirements of statics; namely, 


Is da= Zs, da=P (5) 
Zs = =M (6) 


And in each theory the distribution of strains across the section is assumed 
to be linear. However, in the theory of plasticity it is assumed that the 
stress-strain relationships for the fibers of a member subjected to a simul- 
taneous axial and flexural load may be non-linear and are the same as those 
obtained from simple tension and compression tests. In Figure 2, the stress- 
strain relationship of the material is represented by curve AOB and is as- 
sumed to be the same for each longitudinal fiber. 

The stress in a fiber as given by the theory of plasticity may be interpreted 
as equal to the stress given by the elementary theory plus a correction stress 
to account for the non-linear stress-strain characteristics of the material. 
Thus, 


S=S +S (7) 


The results of the theory of plasticity are compared in this way in Figure 2 to 
those of the elementary theory of mechanics. The ordinates bounded by the 
stress-strain curve AOB are in accordance with the theory of plasticity 
whereas those bounded by the straight line CD are in accordance with the 
elementary theory. The shaded diagram of Figure 2 represents the differ- 
ences in the results of the two theories or correction stresses as defined in 
Equation 7. 

The diagram of correction stresses shown in Figure 2 must satisfy the 
following conditions: 


1. No load or moment is to be contributed by the diagram of correction 
stresses, that is 
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Zs da = 0 (8) 
Zs da = 0 (9) 


2. A linear distribution of strains across the section must be observed. 
3. The correction stresses must account for the non-linear stress-strain 
characteristics of the material. 


In reverse order of above, item 3 is satisfied since the correction 
stresses are represented by the ordinates bounded by the straight line CD 
and the stress-strain curve AOB. Item 2 is satisfied since the abscissas of 
the stress-strain curve and the dimensions of the cross-section are drawn 
to linear horizontal scales. Consequently, the distribution of strains across 
the section is linear. The scale selected for the dimensions of the cross- 
section is in agreement with the values of the strains in the extreme fibers. 
For other values of strains in the extreme fibers than those indicated in the 
Figure, the horizontal scale is either reduced or magnified in proportion to 
the difference in the strains of the extreme fibers. 

The correction stresses indicated in Figure 2 must balance or vanish. 
The terms balance or vanish are used here in the sense that Equations 8 and 


9 are satisfied. The requirements of these Equations are restated as follows: 


1. The sum of the positive correction forces must equal the sum of the nega- 
tive correction forces. A correction force is equal to a correction stress 
times its corresponding differential area. 

2. The centroid of the positive correction forces must coincide with the 
centroid of the negative correction forces. 


The above restatement is particularly useful in pictorial studies of the 
effects of a non-linear stress-strain curve. This should be obvious to those 
already familiar with the pressure line concept for fixed-end beams and 
arches. A procedure of sketching may be developed herein similaz> to the 
procedure used in sketching pressure lines for fixed-end beams.(1) Ina 
fixed-end beam, the angle-changes along the beam caused by lateral loads 
must balance. Thus, conditions of statics for correction forces in a short 
element of a member as discussed herein are analogous to the conditions of 
geometry for angle-changes along a fixed-end beam. If desired, these analo- 
gous relationships may be summarized as in the Column Analogy(2) except 
that the column now has non-linear stress-strain characteristics. In addi- 
tion, the interest herein is in the correction stresses of a column and not in 
the correction moments of a statically indeterminate beam or arch. 

Before pr oceeding to examples of the above procedure, we return to an 
inspection of Figure 2. Two modified moduli of elasticity are introduced; 
namely, 


Ba = (10) 


= As, (11) 


where Sg is the average stress, €g is the strain of the fibers along the cen- 
troidal axis, and Asp and A€ are defined as previously. The values of Eq' 
and Ep’ are reprcsented in the Figure as slopes of corresponding straight 
lines and may de interpreted as reduced moduli of elasticity. They may also 
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be interpreted respectively as measure of stiffness to relative translation 
and to relative rotation of adjacent cross-sections. 

Three examples of the pictorial procedure are given in Figure 3 and are 
discussed below. 

Example A: Given a stress-strain curve, a cross-section, and the values 
of the strains in the extreme fibers. Determine the loads, stresses, angle- 
change per unit of length, and the reduced moduli of elasticity. 

Between the values of the extreme fiber strains, draw a straight line which 
you consider to be the best “fit” to the stress-strain curve. Check “fit” by 
inspecting whether the diagram of correction forces is balanced. If neces- 
sary, revise the straight line until you are satisfied the diagram of correc- 
tion forces is sufficiently balanced. Read in the Figure the values of sq, A4S8p, 
€,, and A€, Now compute the values of load, moment, the reduced moduli of 
elasticity, and the angle-change per unit of length. Stresses in accordance 
with each theory and their differences are observed in the sketch. 

Example B: Given the loads, a cross-section, and the stress-strain curve. 
Determine the stresses, strains, angle-change per unit of length, and the re- 
duced moduli of elasticity. 

Compute the average stress and the stresses in the extreme fibers by 
means of the elementary theory of mechanics. Draw horizontal lines in the 
stress-strain diagram consistent with the computed stresses in the extreme 
fibers. Between these lines, draw a straight line which you consider to be 
the best “fit” to the stress-strain curve. From the intersections of these 
lines drop verticals to the stress-strain curve. The ordinates between the 
straight line and the stress-strain curve represent the correction stresses 
to the elementary theory of mechanics. Check “fit” and, if necessary, revise 
the straight line. Read from the diagrams the value of €,, and 4¢e. Now 
compute the values of the angle-change per unit of length and of the reduced 
moduli of elasticity. 

Example C: Given a stress-strain curve and a cross-section of a pin- 
ended column subjected to a critical axial load. Consider a history of defor- 
mations as in the double- modulus theory(3) wherein an increase in strain of 
a fiber is associated with the tangent modulus and a decrease in strain is 
associated with the initial modulus. Determine the reduced modulus of 
elasticity and the value of L/r which is in agreement with the critical load. 
For illustration only, the bending is assumed to occur about the centroidal 
axis shown in the Figure. 

Compute the average stress on the cross-section and for this value of 
stress draw a horizontal line intersecting the given stress-strain curve. At 
this intersection draw a line tangent to the curve and a line parallel to the 
initial slope of the curve. Consider these two lines as the stress-strain 
curve associated with subsequent buckling. For the critical axial load and 
an arbitrarily assumed moment, compute the stresses in the extreme fibers 
by means of the elementary theory of mechanics. Draw horizontal lines con- 
sistent with the computed stresses in the extreme fibers. Between these 
lines draw a straight line which you consider to be the best “fit” to the 
stress-strain curve associated with buckling. From the intersections of the 
straight line with the two horizontal lines, drop verticals to the assumed 
stress-strain curve. Check “fit” and proceed as in Example B to determine 
the desired position of the straight line. The slope of the straight line is the 
value of the reduced modulus of elasticity. The desired value of L/r is de- 
fined by the modified Euler equation for critical load. 
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For brevity, the above discussion was restricted to that of a cross-section 
having an axis of symmetry and subjected to an axial load and moment about 
an axis normal to the axis of symmetry. The pictorial procedure may be ex- 
tended to the general case of an unsymmetrical cross-section subjected to an 
axial load and to moments about the two principal axes of the section. It may 
be extended to certain problems of torsion, of residual strains, and of other 
histories of deformations than those considered above. It may also be ex- 
tended to studies of the effects of non-uniform changes in temperature or of 
moisture content through the thickness of a long structural member and of 
non-homogeneity through the thickness of the member. Many problems in the 
theory of elasticity which include the effects of warping of cross-sections or 
of non-linear distributions of strains across the sections may be studied in 
the same way. 

If desired, the procedure may be stated more formally and a numerical or 
algebraic procedure may be developed for similar cases as above. 


CONCLUSIONS 


The concept of the pressure line has been shown to be useful for problems 
of an arch with loads in the plane of the arch(1) and for structures curved in 
space. (4) The concept is equally useful for problems involving the inelastic 
deformations of structural elements. 
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Correction Stresses to the 
Elementary Theory 
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| 2 S_ da= 0 
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FIG.2-THE ELEMENTARY THEORY OF MECHANICS 
COMPARED TO A THEORY OF PLASTICITY 
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THE INFINITE ELASTIC BEAM ON A LINEAR 
VISCOELASTIC FOUNDATION 


Alfred M. Freudenthal,! M. ASCE, and Harold G. Lorsch,2 A.M. ASCE 
(Proc. Paper 1158) 


SYNOPSIS 


The infinite, elastic beam on three different types of linear viscoelastic 
foundations is analyzed under various loads. Initial, transient, and terminal 
states of beam deflection and foundation pressure are determined; conditions 
on the load function are stated. Numerical examples with graphs are given. 
The extension of the results to finite beams is discussed. 


Notation 


. coefficients of Fourier series 


p 
. 


. - modulus of elasticity of beam 


- moment of inertia of beam 


. kernel function 


. parameter having the diménsion of length 
. influence length, length of half-period of load 
. integer 


- load per unit length 
. foundation pressure per unit length 


eT . abscissa 


Note: Discussion open until June 1, 1957. Paper 1158 is part of the copyrighted 
Journal of the Engineering Mechanics Division of the American Society of Civil 
Engineers, Vol. 83, No. EM 1, January, 1957. 

1. Prof. of Civ. Eng., Columbia Univ., New York, 


2. Asst. Prof. of Civ. Eng., The City College of New York, N. Y. 
1158-1 


an 
EM1 § 
; 4 
4 
| 
HON 
‘| 
| 


1158-2 EM1 January, 1957 


. reciprocal spring constants 


reciprocal coefficient of viscosity 


T «seeeee24.. ratio between a and ¢, retardation time, relaxation time 


. - coordinate 


Introduction 


History 


The problem of the elastic beam on a continuous foundation has been 
treated extensively in scientific and engineering literature ever since 
Zimmerman(1) gave his solution in 1888. This solution is based on 
Winkler’s(2) assumption that the deflection at any point is proportional to the 
foundation pressure at that point, and does not depend on the pressure at any 
other point of the foundation. 

Wieghart(3) was the first one to investigate the problem under the more 
general assumption that the deflection at any point depends on the foundation 
pressures along a certain length 2L of the foundation. 


+L 


y(x) =const. Jace ak 


(A) 


Here y is the downward deflection of the beam, q the foundation pressure, 
and K a kernel function dependent on the type of foundation. Among others 
Prager,(4) Nemenyi,(5) Marguerre, 6) Biot,(7) Reissner, (8,9) and Volterra(10) 
developed solutions for different kernels K, mostly through the use of Fourier 
integrals. Hetenyi(11) showed, however, that some of these more “accurate” 
solutions in two dimensions may, in effect, be poorer approximations of the 
actual three-dimensional problem than solutions based on the simple Winkler- 
Zimmermann assumption of proportionality between deflection and foundation 
pressure at each point. 

With the recent developments in the field of soil mechanics it became 
evident that, for most types of soils, no kernel function K(| x- |) of the space 
coordinates alone could reproduce the actual behavior of the foundation. It | 
appears necessary to include a function of the time in the relations between 
deflection and foundation pressure. Instead of being considered elastic, the 
foundation should therefore be approximated by a viscoelastic medium; i.e., 
a medium whose force-deformation relations are functions of the time. 

The problem of the elastic beam supported by such a medium may be sim- 
plified by replacing the viscoelastic continuunt with a series of independent, 
viscoelastic elements. This approximation is analogous to the usual engi- 
neering analysis of the elastic beam on an elastic foundation in which the lat- 
ter is replaced by a series of independent, elastic springs. (1,2,11) 


Viscoelastic Models 


The two elements whose combination represents linear viscoelastic be- 
havior are the linearly elastic spring and the viscous dashpot filled with a 
liquid which obeys Newton’s linear law of viscosity.(12) The motion of the 
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perforated piston inside the dashpot produces a resisting force in the liquid 
which is proportional to the velocity of the piston. Figs. 1, 2 and 3 show 
combinations of these elements to form the simplest linear viscoelastic 
models. 

If spring and dashpot are arranged in parallel (Kelvin model, Fig. 1) the 
total force is the sum of the elastic and viscous forces. 


(LA) y4(1/9) (ay/at) (1) 


Here q is the force, y tue deformation, 1/a2 the spring constant, and 1/¢ the 
coefficient of viscosity. Under a constant force q = qo the deformation in- 
creases asymptotically from zero to its final, elastic value. 


(2) 


a-e**) 


where 


t= (3 ) 


is the “retardation time;” it is defined as the time it takes for the deforma- 
tion to reach a-4) of its final, elastic value (Fig. 4). 7 is a measure of the 


delay, caused by the dashpot, in attaining the elastic deformation. On re- 
moval of the force the deformation is gradually recovered. The Kelvin model 
thus represents delayed elasticity or after-effect. 

If spring and dashpot are arranged in series (Maxwell model, Fig. 2), the 
total deformation rate is the sum of the elastic and viscous deformation rates: 


(dy/dt)= (4) 


Under a constant deformation y = y, the initial force Yo/a4 relaxes exponen- 
tially. 


6) 


is the “relaxation time;” it is defined as the time it takes for the force to 
diminish to 1/e of its initial value (Fig. 5). For a constant force the defor- 
mation increases linearly with time. This is the simplest type of the 
phenomenon of creep; thus the Maxwell model represents creep and relaxa- 
tion. 

~~ Fig. 3 shows the arrangement in parallel of a spring and a Maxwell ele- 
ment, frequently referred to as the “Standard Solid.” This model combines 
the two previously mentioned models. If the left spring is very stiff, this 
model tends towards the Kelvin model; if the right spring is very soft the 
model tends towards the Maxwell model. Hence, the Standard Solid model 
exhibits a behavior in between that of a Kelvin and a Maxwell model. 
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Figure 4. 
Deformation-time diagram for a Kelvin 
elewent unaer a constant force. 


Figure 5. 
Force-time diagram for a Maxwell element 
unaer constant deformation (relaxation diagram). 


Solid > 


Figure 6. 
Deformation-time diagrams for viscoelastic 
elements and soil under constant force. 
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Application to Beam on Foundation 


The deformation of a Kelvin element under constant force (Fig. 4) approxi- 
mately corresponds to the consolidation under constant load of a soil, the 
voids of which are filled with water. The behavior of a Maxwell element, on 
the other hand, is somewhat similar to that of an incompressible, cohesive 
soil under shear. The total deformation of a cohesive soil is a combination 
of deformations due to consolidation and shear. Consolidation is the pre- 
dominant long-term effect under relatively light loads, while shear becomes 
more significant with increasing load intensities, particularly over short time 
intervals. Fig. 6 shows how the Standard Solid modei best approximates ac- 
tual soil behavior, while the Kelvin foundation approximates long-term be- 
havior and the Maxwell foundation short-term behavior. 

The present paper treats an elastic beam on various viscoelastic founda- 
tions corresponding to models of Figs. 1, 2, and 3. In all cases it will be as- 
sumed that the foundation pressure and its rate of change with time at a 
point depend on the deflection and its rate of change with time at that point 
only; in other words, the foundation will be idealized as a series of identical, 
viscoelastic elements, spaced uniformly under the beam. 


Elastic Beam on a Linear Viscoelastic (Standard Solid) Foundation 
With the notation of Fig. 3 it follows from (4) and from y =ao9q9, using 
= qj + Gg, that 
This becomes 
pt) (dq/at), (7) 
if the notation 


is used. The solution of this differential equation is 


t 
y= (4192/9) fae“ at) (9) 


which gives the deformation y as a function of a given force q. 

The differential equation of the elastic beam on a viscoelastic foundation, 
loaded by an arbitrary load p(x,t) (Fig. 7), is obtained by substitution of (9) 
into the differential equation of the elastic beam(13) 


B1 )ap (x, t)-q (x,t) (10) 


where E and I are modulus of elasticity and moment of inertia, respectively, 
of the beam. By differentiation with respect to t and subsequent elimination 
of the integral in (9), the following equation is obtained. 
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p(x) 


Figure 7. 
Elastic beam on a linear viscoelastic 
(Standard Solid) foundation. 


P(X) 


Figure 8, 
Infinite elastic beam on a foundation 
with after-effect (Kelvin foundation). 
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a +k +2848. 2,28 
atax Ox* at (12) 
Here 
Qytag) = Ela, a o/pt (12) 


where k is a length. Whenever (11) is solved for the soil pressure q(x,t), the 
deflection of the beam y(x,t) can be obtained from (9). A particularly simple 
solution arises in the case of (84p/a x4) = 0 for which 


a(x, tmp (x,t)+ (A, cosh sinh mx 
meO 


+ C,cos mx+D sin mx) e + I/mtk# 


Further, simpler solutions will now be given for the case of the elastic beam 
on either the Kelvin or the Maxwell foundation. 


Elastic Beam on a Foundation with After-Effect (Kelvin Foundation) 


General Formulae 


If the left spring in the model of Fig. 3 is made infinitely stiff (a; = 0), the 
model changes into that of Fig. 1, and the beam of Fig. 7 becomes the beam 
on a Kelvin foundation (Fig. 8). Its differential equation is obtained by sub- 
stituting (1) into (10). 


34 Ov(x,t 
« 

Here 7 is defined by (3) and k by 
4 
k = Ela,. 

The initial condition for this beam is that the deflection is zero att = 0. The 
deflection of a beam on an elastic foundation is a particular solution of (14). 
The solution of the homogeneous equation for (14) (right-hand side equal to 
zero) is obtained by setting y(x,t) = f(x) - g(t) which splits (14) into two ordi- 
nary differential equations that can readily be solved. 


4,4 
(x,t) = (Agcosh mx¢B,sinh mx+C,cos m+D,sin mx) (16) 


where m is a non negative integer. Hence the general solution of (14) is the 
sum of solutions (16) and the deflection of the beam on an elastic foundation. 
The soil pressure is obtained by substituting this solution in (1). 
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Infinite Beam Loaded by a Space-periodic Load, Invariant with Time 


If the beam is loaded by a space-periodic load of period 27* which does 
not vary with time (Fig. 9a), the load function can be expanded in a Fourier 
series 


mae 


P(x) af/2 + (a, cos mx+b, sin mx) (17) 


subject to very general conditions.(14) The coefficients of this series are 
given by 


a, = cm fre cos mx dx (17a) 
by = (1/m) fre sin mx dx (17b) 


The corresponding elastic foundation pressures and deflections are 


maco 


Yel. @ 8/2 + & (1/Litmt«t J ) (ag cos mx+b, sin mx) (18a) 


Yel. (18b) 


as can be verified from (10). Eqs. (18) thus furnish the particular solution 
for (14). The solution of the homogeneous equation is given by (16) with the 
following substitutions. 


Q, Cowan /2 By" 0, 
44 4,4) 1.2.2 


Hence, the complete solution for the deflection is 


y (x, t)*(Qpa0/2) (1-07 COs MX (49) 


+b, sin mJ[l-e~ 


The foundation pressure follows by substitution in (1). 


mec 


q(x, +> (1/C1+m*«*I) (a, cos mx (20) 
me| 


4.4 
sin mx) (1+m*«*e~ (1+m*« 
* If the period of the load is 2L instead of 27, m is to be replaced by mz/L 


in all the equations that follow, except for the limits of all the sums which 
remain unchanged from m = 1 to m =oo, 
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Eqs. (19) and (20) clearly show the behavior of the beam. At the instant of 
application of the load (t = 0) the foundation pressure is equal to the load, and 
the deflection is zero; the beam remains stiff at this instant, and the load is 
transmitted vertically through the beam and counteracted by an equal and 
opposite foundation pressure: q(x,0) = p(x,0).* 

After an infinite time interval (t>~) foundation pressure and deflection ap- 
proach their elastic values. During the transient state (0St <= ) the 
approach to the final state for the beam on the Kelvin foundation is the faster, 
the larger m, k, and 1/7; i.e. the smaller the period of the load (see footnote, 
page 8), the stiffer the beam with respect to the foundation, and the lower the 
viscosity of the foundation. ; 


Infinite Beam Loaded by an Arbitrary Load, Invariant with Time 


If the infinite, elastic beam is loaded by an arbitrary load which does not 
vary with time, the load function can be represented by its Fourier integral. 


+= 
p(x) =(1/m) f amfp(z) cos [m(z-x)] dz (21) 


subject to very general conditions.(15) By direct analogy with the procedures ’ 
used in the previous paragraph the elastic foundation pressures and deflec- 
tions are 


+ 
gers = f cos{m(z-x)] dz (22a) 
—@ 


Gy) anf (1+m J p(z) cos[m(z-xJaz (22) 


The last equation again furnishes a particular solution for (14). The solution 
of the homogeneous equation is taken as the integral over m of solutions(16) 
setting Am=Bm=Dm=0 and determining Cy, such that the deflection vanishes 
at t=0. 


y(x,t)= (afm) f (1-e~ te) ony (14m (23) 

fp(z) cos{m(z-x)] az. 

As before, the foundation pressure follows from (1). i 

q(x, t)=(1/m) f (14m (14m amy 
° 
(24) 


cos[m(z-x)] dz. 


* This is true, provided the spreading out of the load through the thickness of 
the beam is not considered. Such effects are taken into account in 
reference 19. 
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Cb) 


(Cc) 


Figure 9, 
a) Beam loaded by a 4 oe a load p(x); 
bd) Average loaa Py a 
c) Self-equiiibrating ldad 


Figure 10, 
Infinite beam loaded by a uniformly 
aistributed load over a finite region. 
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The interpretation of the beam behavior is similar to that discussed above 
for the space-periodic load. 


Infinite Beam Loaded by a Uniformly Distributed Load over a Finite Region 


Let the infinite, elastic beam of Fig. 10 be loaded by a uniformly distribut- 
ed, time invariant load pp over the finite region -b&xS+b. Then the second 
integrals in eqs. (21) through (24) extend over the limits -b§$z§+b only, and 
the integration can be readily performed: 


<A+n*k*)tA 
y (x, (ef) | (1-e Jam (25) 


and 
q(x, t)=(2p,/m) [stn gon (26) 
m(1+m"k*) 


After an infinite time interval the parenthesis in (25) becomes equal to unity, 
and 


is the Fourier integral of the expression for the deflection of a beam on an 
elastic foundation, generally found in the literature.(16) As a rule, two equa- 
tions for ye), are given depending on whether the point under consideration is 
within, or outside of, the loaded region. Eq. (27) covers both cases with a 

single expression. 


Infinite Beam Loaded by a Single, Concentrated, Constant Load 


If the infinite, elastic beam is acted upon by a single, concentrated, time 
invariant load, the preceding Fourier integral solutions are applicable with 
the following substitutions: 2bp>P and b>0. Eqs. (25) and (26) then become, 
respectively, 


4,4 
y (x,t) f t/t) am (28) 


and 


4,4 
© 


The last integral does not converge at t = 0 on account of the sharp dis- 
continuity of load and foundation pressure at the origin at that instant. This, 
however, is not a serious defect of the solution because it is known that, at 
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t = 0, the foundation pressure equals the load.* 
For t>~= eq. (28) can be written as 


This is the Fourier integral of 


Yel. = (0, P/2 (x/¥2k)+sin (x, (31) 


which is the deformation of a beam on an elastic foundation, loaded by a 
single concentrated load P.(17) This again shows the general validity of the 
remarks concerning the behavior of an elastic beam on a (Kelvin) foundation 
with after-effect which were made under (V,b) above. 


Elastic Beam on a Relaxing (Maxwell) Foundation 


General Formulae 


If the right spring in the model of Fig. 3 is eliminated (a9), this model 
changes into that of Fig. 2, and the beam of Fig. 7 becomes the beam on a 
Maxwell foundation (Fig. 11). Eq. (11) simplifies to 


Figure ll. 
Infinite elastic beam on a 
relaxing (Maxwell) foundation. 


*Also see footnote, page 10. 
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where T is defined by (6) and k by 


x! = Ela, . (33) 


If the applied loads do not vary with time, (32) reduces to 


49° t x t 


Solutions of this equation are assumed to be of the form q(x,t) = f(x)-g(t) 
which splits (32a) into two ordinary differential equations. 


Qn #(A,cosh mx+B, sinh mx4¢C,cos mx 


+D,sin mx) e { Tat } 


Here m is a non-negative integer. The deflection y(x,t) is obtained from the 
foundation pressure q(x,t) by integration of (4) with respect to the time. 


t 


y (x,t) = at +F(x) 
° 


where F(x) is an arbitrary function of x alone. It can be evaluated by con- 
sidering that, at t = 0, the deflection equals the elastic deflection. 


(x, 0)=F (x)= 4(x,0). 


The use of the last two equations yields the deformation 


4 
y(x,t)=aa (x,t)+ of q(x,t) at 


after q(x,t) has been obtained from (34). 


Infinite Beam Loaded by a Space-periodic Load, Invariant with Time 


If the infinite beam is loaded by a space-periodic load (Fig. 9a) which does 
not vary with time, the Fourier series expansion (17) can again be used. The 
foundation pressures are obtained by an infinite sum of solutions (34) with the 
following substitutions. 


Ao™Bo®Do=0, Co®ag/2 A,=B= 0, 


4 
(34) 
i 
fag (35) 
4 
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Hence the Maxwell foundation pressure equals 
-t/e 


i+m?k* 


and the beam deflection follows from (35). 


m=co 


y (x,t) =(qa,/2) fe) +> (a, /n*k*) (a cos mx 


-t/t } 
1+1/m*k* (37) 
€ 
+ b,sin mx) (1 - ) 

The behavior of the beam becomes apparent from the last equation: at the 
instant of load application (t=0) the deflection equals the elastic deflection 
(18b). After an infinite time interval (t>~) the beam deflection is composed 
of two parts: 

I. The deflection of a beam loaded by a uniform load equal to the average 
load pp, = 49/2 (Fig. 9b); this is the uniform settlement of a straight, rigid 
beam corresponding to the uniform motion of the piston in the viscous dash- 
pot (Fig. 2) under a constant force; 

II.‘ the elastic deflection of the beam under a self-equilibrating force sys- 
tem equal to the difference between the actual load and its average (Fig. 9c). 
This latter deflection is entirely independent of the stiffness and viscosity of 
the foundation, and only depends on the bending stiffness EI of the beam. 

The final state of the beam, therefore, is that of a floating, elastic beam. 
During the transient state (OStS e) the approach to the final deflection is the 
faster the larger m, k, and1/7. The same criteria were found to hold for 
the beam on the Kelvin foundation. From (36) it is apparent that the founda- 
tion pressure at t = 0 represents the elastic pressure (18a) and that, with in- 
finite time, it approaches a constant value equal to the average load po = ao /2. 


Infinite Beam Loaded by an Arbitrary Load, Invariant with Time 


For the infinite, elastic beam loaded by a time invariant, arbitrary load 
the Fourier integral representation of the load (21), elastic foundation pres- 
sure (22a), and elastic deflection (22b) may again be used. The Maxwell foun- 
dation pressure follows bv integration over m of solutions (34) with 
Am=Bm=Dm=0. Cm is determined from the initial (elastic) state. 


a(x,t) = J p(z) cos[m(z-x)Jdz (38) 


1lem*k* 


The integration (35) yields the deformation. 


+00 


y (x, t)=(a, f (1- fe) cos[m(z-x)] az (39) 


i. 
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The variation of deflections and foundation pressures is identical to that 
discussed in the preceding paragraph; i.e., the initial state is that of a beam 
on an elastic foundation, the final state is that of a floating beam. Eqs. (38) 
and (39), however, are valid solutions only if, for any point x, a finite interval 
X1< xX <xg exists such that 


[ p(x)-po J axed (40) 


where po is the average load on the beam. 
If equation (40) is not satisfied, the beam becomes subject to an infinitely 
large bending moment as t approaches infinity. Obviously, no real structure 
can sustain infinite moments, and a redistribution of stresses will take place 
once the proportional limit of the material in the elastic beam has been ex- 
ceedeii. 

Among others, the regional load discussed in (V,d) above and the single 
concentrated load (V,e) do not satisfy (40). 


Numerical Examples 


A reinforced concrete beam, 30 in. wide and 20 in. deep, is continuously 
supported by a Kelvin or a Maxwell foundation. The beam is loaded by 30 x 
30 in. columns every 20 ft. throughout its infinite length (Fig. 12). Each 
column load of 300 kips is considered uniformly distributed over the contact 
area. The modulus of elasticity of the concrete is 2500 kips per sq. in., the 
foundation modulus 0.75 kips per sq. in. per in. The last constant must be 
multiplied by the beam width in order to obtain 1/a, the reciprocal of the 


20'-on 201-0n 


201-0" 


3do* 300 2do* 


Maxwell or Kelv oundation 


Figure 12. 
Elastic beam on a viscoelastic foundation 
loaded by equal column loads. 
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elastic foundation constant used in this paper. This is due to the fact that the 
foundation modulus is defined as the ratio of bearing unit stress to deflection, 
while o@ is the ratio of deflection to bearing load per unit length. 

The Fourier coefficients for the load function p(x) are determined from 
(17a) and (17b), 


m/2 sin(m7/16) .... if m is even 


... if m is odd (41) 


m = (40/mz) (-1) 


The deflections and bearing pressures for the beam on the Kelvin and Max- 
well foundations, respectively, are obtained from equations (19), (20), (36), 
and (37) when m is replaced by m7/L, the summations are extended over all 
even integers, and relations (41) are used. Bending moments in the beam are 
obtained through term by term differentiation of equations (19) and (37). Suf- 
ficient numbers of terms of the series were taken to make the first term 
omitted less than 1% of the sum of the previous terms. This somewhat arbi- 
trary criterion limits the error to approximately 2% for series, all of whose 
terms are positive, and to a little over 1/2% for series with alternating posi- 
tive and negative terms. The former applies to points on the beam under the 
columns, the latter to all other points calculated. As a rule, six terms were 
thus required for the bending moments and three terms for the bearing pres- 
sures and deflections. 

Figs. 13 through 15 show the variations in space and time of the bearing 
pressure, bending moment, and deflection for the beam on a foundation with 
after-effect (Kelvin foundation). The corresponding data for the beam on a 
relaxing (Maxwell) foundation are given in Figs. 16 through 18. The variation 
with time of the beam behavior is shown in terms of the relative time t/T. 
This makes Figs. 13 through 18 valid for any value of ¢, provided the elastic 
properties of beam and foundation are those assumed above. Bearing pres- 
sures and bending moments for the beam on a viscoelastic foundation can be 
expressed in terms of the bearing pressures and bending moments for the 
beam on an elastic foundation, if the viscous property of the foundation (¢ or 
T) and the ratio of the elastic properties of beam and foundation (k) are given. 
The actual values of those properties must be known, however, in order to ob- 
tain values for the deflections. 

Fig. 13 shows the gradual spreading-out of the bearing pressure for the 
beam on the foundation with after-effect (Kelvin foundation) from the state of 
the rigid foundation to that of the elastic foundation. Small negative bearing 
pressures occur in a region between the columns during a short time after the 
load has been applied. Even though these small tensile stresses exceed the 
adhesion between beam and foundation and hence the actual bearing stress 
distribution changes, this change will be limited to a short initial time inter- 
val. The corresponding changes in bending moments and deflections will 
likewise be limited to a short period of time and will not change the long-term 
behavior of the beam. The upward deflections over a small portion of the 
beam for t/r less than 2 (Fig. 15) are, in part, a consequence of those tensile 
stresses. However, even though no tensile stress between beam and founda- 
tion can be transmitted in a particular foundation, these upward deflections 
will not disappear completely. 
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(in kips per inch) 


for beam on Kelvin 


Figure 15. x t = 

Bearing pressure é 
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Figure 15. 


Deflection (in inches) for beam on Kelvin foundation. 
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Fig. 14 indicates that, at some time or other, negative moments occur al- 
most throughout the length of the beam on a Kelvin foundation. These mo- 
ments are small, compared to moments for the beam on an elastic foundation 
(t/7r =~), but they make it advisable to provide a nominal amount of top rein- 
forcing steel over the entire length of a reinforced concrete beam. This is 
generally done in beams where such steel is needed for the proper placement 
of stirrups, but it is frequently omitted in beams for which no stirrups are 
required. The curves in Figs. 13 through 15 are not materially affected by a 
change in relative stiffness between »eam and foundation. 

The gradual flattening-out of the bearing pressure curve for a beam ona 
relaxing (Maxwell) foundation from its elastic to its average value is shown 
in Fig. 16. Fig. 17 shows how the bending moment in the beam increases as 
a consequence. Both bearing pressure and bending moment remain virtually 
constant at the quarter-points. The maximum positive bending moment in- 
creases by 57% over that for the beam on an elastic foundation, the maximum 
negative bending moment by 138%. This shows that, as far as stresses in the 
beam are concerned, it is unsafe to assume elastic behavior for a foundation 
in which irrecoverable shear deformations play a significant part. 

If the elastic stiffness of the foundation with respect to that of the beam is 
decreased (k increases), bearing pressures and bending moments for the 
beam on an elastic foundation gradually approach those for the beam on a re- 
laxing foundation at infinite time; i.e., those of the floating beam. For a value 
of k = 10 feet, their maximum difference is about 2%. (Figs. 13 through 18 
are based on k = 3.22 feet) These differences become larger with increasing 
relative stiffness of the foundation and tend to become infinite for a perfectly 
flexible beam or an infinitely stiff foundation: k = 0 (Fig. 19). 

Fig. 18 illustrates the steady downward motion as a whole of the beam on a 
Maxwell foundation on which an increase in beam flexure is superposed. An 
increase in k (more flexible foundation or stiffer beam) makes all the curves 
flatter, a decrease in k has the opposite result; neither changes the essential 
two-phase character of the deflection curves. 


SUMMARY 


A general differential equation is derived for an elastic beam on a linear 
viscoelastic foundation, acted upon by arbitrarily varying loads in space and 
time (11). A solution is given for the special case of (04p/a x4) = 0. Eq. (13). 

Assuming less general viscoelastic laws for the foundation, the time and 
space variations of foundation pressure and deflection are determined for an 
infinite elastic beam on a foundation with linear after-effect (Kelvin founda- 
tion)(1) and a linear relaxing (Maxwell) foundation(4) provided the loads do 
not vary with time. The initial state of the beam on a Kelvin foundation is 
the beam on a perfectly rigid foundation, the final state of the beam ona 
Kelvin foundation and the initial state of the beam on a Maxwell foundation is 
the beam on an elastic foundation, while the final state of the beam on a Max- 
well foundation is the floating beam. 

For the infinite beam on a relaxing (Maxwell) foundation the load function 
must satisfy condition (40); otherwise the moments in the beam will become 
infinite if the assumption of linear behavior for the beam material is main- 
tained. Among others, the single concentrated load does not satisfy equation 
(40). 
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Figure i6. 


Bearing pressure (in kips per inch) 


for beam on Maxwell foundation, 
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Figure 18. 
Deflection (in inches) for beam on Maxwell foundation. 
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Figure 19. 
Ratio of maxinun negative bending moment in beam on Maxwell 
foundation to that in a beam on an elastic foundation as 
a function of relative stiffness of beam and foundation, 


ASCE FREUDENTHAL - LORSCH 1158-21 
column 
0.2 
t 0.2 
= 
\ 
ight 
Z 
#=20 
1.8 
elastic 
| | 
| 
| 
5 5 10 


1158-22 EM1 January, 1957 


Graphs are given which show the behavior of an infinite, elastic beam un- 
der equally spaced column loads on a foundation with after-effect and a 
relaxing foundation. 

All the differential equations treated in this paper are linear; hence, any 
number of solutions for a particular foundation may be used in superposition 
to obtain additional solutions. The usual methods of transition from infinite 
to finite beams by the use of end forces and.end moments(18) is also valid 
for beams on a viscoelastic foundation. In general, however, these end 
forces and end moments are functions of the time. 
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